home *** CD-ROM | disk | FTP | other *** search
/ Languguage OS 2 / Languguage OS II Version 10-94 (Knowledge Media)(1994).ISO / language / asa / asa.c < prev    next >
C/C++ Source or Header  |  1994-10-23  |  93KB  |  3,163 lines
  1. /***********************************************************************
  2. * Lester Ingber (copyright) (c)
  3. * See COPYING License in this directory
  4. * Date   1 Jan 93
  5. ***********************************************************************/
  6.  
  7. #define ASA_ID \
  8. "/* $Id: asa.c,v 4.2 1994/10/23 23:35:16 ingber Exp ingber $ */"
  9.  
  10. #include "asa.h"
  11.  
  12. /* static variables in EXPONENT_CHECK () */
  13. static double MIN_EXPONENT, MAX_EXPONENT;
  14.  
  15. /***********************************************************************
  16. * asa
  17. *       This procedure implements the full ASA function optimization.
  18. ***********************************************************************/
  19. #if HAVE_ANSI
  20. double
  21. asa (double (*user_cost_function) (
  22.                double *, double *, double *, double *, double *,
  23.               ALLOC_INT *, int *, int *, int *, USER_DEFINES *),
  24.      double (*user_random_generator) (),
  25.      double *parameter_initial_final,
  26.      double *parameter_minimum,
  27.      double *parameter_maximum,
  28.      double *tangents,
  29.      double *curvature,
  30.      ALLOC_INT * number_parameters,
  31.      int *parameter_type,
  32.      int *valid_state_generated_flag,
  33.      int *exit_status,
  34.      USER_DEFINES * OPTIONS)
  35. #else
  36. double
  37. asa (user_cost_function,
  38.      user_random_generator,
  39.      parameter_initial_final,
  40.      parameter_minimum,
  41.      parameter_maximum,
  42.      tangents,
  43.      curvature,
  44.      number_parameters,
  45.      parameter_type,
  46.      valid_state_generated_flag,
  47.      exit_status,
  48.      OPTIONS)
  49.      double (*user_cost_function) ();
  50.      double (*user_random_generator) ();
  51.      double *parameter_initial_final;
  52.      double *parameter_minimum;
  53.      double *parameter_maximum;
  54.      double *tangents;
  55.      double *curvature;
  56.      ALLOC_INT *number_parameters;
  57.      int *parameter_type;
  58.      int *valid_state_generated_flag;
  59.      int *exit_status;
  60.      USER_DEFINES *OPTIONS;
  61. #endif
  62. {
  63.   int index_cost_constraint,    /* index initial cost functions averaged
  64.                    to set initial temperature */
  65.     tmp_var_int,        /* integer temporary value */
  66.     index_cost_repeat;        /* test OPTIONS->COST_PRECISION when =
  67.                    OPTIONS->MAXIMUM_COST_REPEAT */
  68.  
  69.   ALLOC_INT index_v,        /* iteration index */
  70.    *start_sequence;        /* initial OPTIONS->SEQUENTIAL_PARAMETERS
  71.                    used if >= 0 */
  72.   double sampled_cost_sum,    /* temporary initial costs */
  73.     final_cost,            /* best cost to return to user */
  74.     tmp_var_db1, tmp_var_db2;    /* temporary doubles */
  75.   int *curvature_flag;
  76.   FILE *ptr_asa_out;        /* file ptr to output file */
  77.  
  78.   /* The 3 states that are kept track of during the annealing process */
  79.   STATE *current_generated_state, *last_saved_state, *best_generated_state;
  80.  
  81. #if ASA_PARALLEL
  82.   /* parallel generated states */
  83.   STATE *gener_block_state;
  84. #endif
  85.  
  86.   /* Holds values of parameter increments */
  87.   double *delta_parameter;
  88.   /* Holds values of quench_param scales */
  89.   double *quench_param;
  90.   /* Holds value of quench_cost scales */
  91.   double *quench_cost;
  92.  
  93.   /* The array of tangents (absolute value of the numerical derivatives),
  94.      and the maximum |tangent| of the array */
  95.   double *maximum_tangent;
  96.  
  97.   /* ratio of acceptances to generated points - determines when to
  98.      test/reanneal */
  99.   double *accepted_to_generated_ratio;
  100.  
  101.   /* temperature parameters */
  102.   double temperature_scale, *temperature_scale_parameters;
  103.   /* relative scalings of cost and parameters to temperature_scale */
  104.   double *temperature_scale_cost;
  105.   double *current_user_parameter_temp;
  106.   double *initial_user_parameter_temp;
  107.   double *current_cost_temperature;
  108.   double *initial_cost_temperature;
  109.   double log_new_temperature_ratio;    /* current *temp = initial *temp *
  110.                        exp(log_new_temperature_ratio) */
  111.   ALLOC_INT *index_exit_v;    /* information for asa_exit */
  112.  
  113.   /* counts of generated states and acceptances */
  114.   LONG_INT *index_parameter_generations;
  115.   LONG_INT *number_generated, *best_number_generated_saved;
  116.   LONG_INT *recent_number_generated, *number_accepted;
  117.   LONG_INT *recent_number_acceptances, *index_cost_acceptances;
  118.   LONG_INT *number_acceptances_saved, *best_number_accepted_saved;
  119.  
  120.   /* Flag indicates that the parameters generated were
  121.      invalid according to the cost function validity criteria. */
  122.   LONG_INT *number_invalid_generated_states;
  123.   LONG_INT repeated_invalid_states;
  124.  
  125. #if ASA_PARALLEL
  126.   LONG_INT index_parallel;    /* count of parallel generated states */
  127.   LONG_INT parallel_generated;    /* saved *recent_number_generated */
  128.   LONG_INT parallel_block_max;    /* saved OPTIONS->gener_block_max */
  129. #endif
  130.  
  131.   /* used to index repeated and recursive calls to asa */
  132.   /* This assumes that multiple calls (>= 1) _or_ recursive
  133.      calls are being made to asa */
  134.   static int asa_open = FALSE;
  135.   static int number_asa_open = 0;
  136.   static int recursive_asa_open = 0;
  137.  
  138.   /* static variables used in EXPONENT_CHECK() */
  139.   MIN_EXPONENT = 0.9 * log ((double) MIN_DOUBLE);
  140.   MAX_EXPONENT = 0.9 * log ((double) MAX_DOUBLE);
  141.  
  142.   /* initializations */
  143.  
  144.   if ((curvature_flag =
  145.        (int *) calloc (1, sizeof (int))) == NULL)
  146.       exit (9);
  147.   if ((maximum_tangent =
  148.        (double *) calloc (1, sizeof (double))) == NULL)
  149.       exit (9);
  150.   if ((accepted_to_generated_ratio =
  151.        (double *) calloc (1, sizeof (double))) == NULL)
  152.       exit (9);
  153.   if ((temperature_scale_cost =
  154.        (double *) calloc (1, sizeof (double))) == NULL)
  155.       exit (9);
  156.   if ((current_cost_temperature =
  157.        (double *) calloc (1, sizeof (double))) == NULL)
  158.       exit (9);
  159.   if ((initial_cost_temperature =
  160.        (double *) calloc (1, sizeof (double))) == NULL)
  161.       exit (9);
  162.   if ((index_exit_v =
  163.        (ALLOC_INT *) calloc (1, sizeof (ALLOC_INT))) == NULL)
  164.     exit (9);
  165.   if ((start_sequence =
  166.        (ALLOC_INT *) calloc (1, sizeof (ALLOC_INT))) == NULL)
  167.     exit (9);
  168.   if ((number_generated =
  169.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  170.     exit (9);
  171.   if ((best_number_generated_saved =
  172.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  173.     exit (9);
  174.   if ((recent_number_generated =
  175.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  176.     exit (9);
  177.   if ((number_accepted =
  178.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  179.     exit (9);
  180.   if ((recent_number_acceptances =
  181.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  182.     exit (9);
  183.   if ((index_cost_acceptances =
  184.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  185.     exit (9);
  186.   if ((number_acceptances_saved =
  187.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  188.     exit (9);
  189.   if ((best_number_accepted_saved =
  190.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  191.     exit (9);
  192.   if ((number_invalid_generated_states =
  193.        (LONG_INT *) calloc (1, sizeof (LONG_INT))) == NULL)
  194.     exit (9);
  195.  
  196.   if ((current_generated_state =
  197.        (STATE *) calloc (1, sizeof (STATE))) == NULL)
  198.     exit (9);
  199.   if ((last_saved_state =
  200.        (STATE *) calloc (1, sizeof (STATE))) == NULL)
  201.     exit (9);
  202.   if ((best_generated_state =
  203.        (STATE *) calloc (1, sizeof (STATE))) == NULL)
  204.     exit (9);
  205. #if ASA_PARALLEL
  206.   if ((gener_block_state =
  207.        (STATE *) calloc (OPTIONS->gener_block_max, sizeof (STATE))) == NULL)
  208.     exit (9);
  209. #endif
  210.  
  211.   if (asa_open == FALSE)
  212.     {
  213.       asa_open = TRUE;
  214.       ++number_asa_open;
  215. #if ASA_PRINT
  216.       if (number_asa_open == 1)
  217.     {
  218.       /* open the output file */
  219. #if USER_ASA_OUT
  220.       if (!strcmp (OPTIONS->asa_out_file, "STDOUT"))
  221.         {
  222.           ptr_asa_out = stdout;
  223.         }
  224.       else
  225.         {
  226.           ptr_asa_out = fopen (OPTIONS->asa_out_file, "w");
  227.         }
  228. #else
  229.       if (!strcmp (ASA_OUT, "STDOUT"))
  230.         {
  231.           ptr_asa_out = stdout;
  232.         }
  233.       else
  234.         {
  235.           ptr_asa_out = fopen (ASA_OUT, "w");
  236.         }
  237. #endif
  238.     }
  239.       else
  240.     {
  241. #if USER_ASA_OUT
  242.       ptr_asa_out = fopen (OPTIONS->asa_out_file, "a");
  243. #else
  244.       ptr_asa_out = fopen (ASA_OUT, "a");
  245. #endif
  246.       fprintf (ptr_asa_out, "\n\n\t\t number_asa_open = %d\n",
  247.            number_asa_open);
  248.     }
  249. #endif
  250.     }
  251.   else
  252.     {
  253.       ++recursive_asa_open;
  254. #if ASA_PRINT
  255.       if (recursive_asa_open == 1)
  256.     {
  257.       /* open the output file */
  258. #if USER_ASA_OUT
  259.       ptr_asa_out = fopen (OPTIONS->asa_out_file, "w");
  260. #else
  261.       ptr_asa_out = fopen (ASA_OUT, "w");
  262. #endif
  263.     }
  264.       else
  265.     {
  266.  
  267. #if USER_ASA_OUT
  268.       ptr_asa_out = fopen (OPTIONS->asa_out_file, "a");
  269. #else
  270.       ptr_asa_out = fopen (ASA_OUT, "a");
  271. #endif
  272.       fprintf (ptr_asa_out, "\n\n\t\t recursive_asa_open = %d\n",
  273.            recursive_asa_open);
  274.     }
  275. #endif
  276.     }
  277.  
  278. #if ASA_PRINT
  279.   /* print header information as defined by user */
  280.   print_asa_options (ptr_asa_out, OPTIONS);
  281.  
  282. #if TIME_CALC
  283. /* print starting time */
  284.   print_time ("start_asa", ptr_asa_out);
  285. #endif
  286.   fflush (ptr_asa_out);
  287. #endif
  288.  
  289. /* set indices and counts to 0 */
  290.   *best_number_generated_saved =
  291.     *number_generated =
  292.     *recent_number_generated =
  293.     *recent_number_acceptances = 0;
  294.   *index_cost_acceptances =
  295.     *best_number_accepted_saved =
  296.     *number_accepted =
  297.     *number_acceptances_saved = 0;
  298.   index_cost_repeat = 0;
  299.  
  300. #if ASA_SAMPLE
  301.   OPTIONS->n_accepted = 0;
  302.   OPTIONS->average_weights = 1.0;
  303. #endif
  304.  
  305. /* do not calculate curvatures initially */
  306.   *curvature_flag = FALSE;
  307.  
  308. /* allocate storage for all parameters */
  309.   if ((current_generated_state->parameter =
  310.        (double *) calloc (*number_parameters, sizeof (double))
  311.       ) == NULL)
  312.       exit (9);
  313.   if ((last_saved_state->parameter =
  314.        (double *) calloc (*number_parameters, sizeof (double))
  315.       ) == NULL)
  316.       exit (9);
  317.   if ((best_generated_state->parameter =
  318.        (double *) calloc (*number_parameters, sizeof (double))
  319.       ) == NULL)
  320.       exit (9);
  321. #if ASA_PARALLEL
  322.   parallel_block_max = OPTIONS->gener_block_max;
  323.   parallel_generated = OPTIONS->gener_block;
  324.  
  325.   for (index_parallel = 0; index_parallel < parallel_block_max;
  326.        ++index_parallel)
  327.     {
  328.       if ((gener_block_state[index_parallel].parameter =
  329.        (double *) calloc (*number_parameters, sizeof (double))
  330.       ) == NULL)
  331.       exit (9);
  332.     }
  333. #endif
  334.  
  335.   if ((initial_user_parameter_temp =
  336.        (double *) calloc (*number_parameters, sizeof (double))
  337.       ) == NULL)
  338.       exit (9);
  339.   if ((index_parameter_generations =
  340.        (LONG_INT *) calloc (*number_parameters, sizeof (LONG_INT))
  341.       ) == NULL)
  342.     exit (9);
  343.  
  344. /* set all delta_parameter[] */
  345.   if (OPTIONS->DELTA_PARAMETERS == TRUE)
  346.     {
  347.       delta_parameter = OPTIONS->user_delta_parameter;
  348.       VFOR (index_v)
  349.     delta_parameter[index_v] =
  350.     OPTIONS->user_delta_parameter[index_v];
  351.     }
  352.   else
  353.     {
  354.       if ((delta_parameter =
  355.        (double *) calloc (*number_parameters, sizeof (double))
  356.       ) == NULL)
  357.       exit (9);
  358.       VFOR (index_v)
  359.     delta_parameter[index_v] = OPTIONS->DELTA_X;
  360.     }
  361.  
  362. /* set all quench_param[] */
  363.   if (OPTIONS->QUENCH_PARAMETERS == TRUE)
  364.     {
  365.       quench_param = OPTIONS->user_quench_param_scale;
  366.       VFOR (index_v)
  367.     quench_param[index_v] =
  368.     OPTIONS->user_quench_param_scale[index_v];
  369.     }
  370.   else
  371.     {
  372.       if ((quench_param =
  373.        (double *) calloc (*number_parameters, sizeof (double))
  374.       ) == NULL)
  375.       exit (9);
  376.       VFOR (index_v)
  377.     quench_param[index_v] = ONE;
  378.     }
  379.  
  380. /* set quench_cost[] */
  381.   if (OPTIONS->QUENCH_COST == TRUE)
  382.     {
  383.       quench_cost = OPTIONS->user_quench_cost_scale;
  384.       *quench_cost = OPTIONS->user_quench_cost_scale[0];
  385.     }
  386.   else
  387.     {
  388.       if ((quench_cost =
  389.        (double *) calloc (1, sizeof (double))) == NULL)
  390.       exit (9);
  391.       *quench_cost = ONE;
  392.     }
  393.  
  394. /* set all temperatures */
  395.   if ((current_user_parameter_temp =
  396.        (double *) calloc (*number_parameters, sizeof (double))
  397.       ) == NULL)
  398.       exit (9);
  399.   if (OPTIONS->USER_INITIAL_PARAMETERS_TEMPS == TRUE)
  400.     {
  401.       VFOR (index_v)
  402.     current_user_parameter_temp[index_v] =
  403.     initial_user_parameter_temp[index_v] =
  404.     OPTIONS->user_parameter_temperature[index_v];
  405.     }
  406.   else
  407.     {
  408.       VFOR (index_v)
  409.     current_user_parameter_temp[index_v] =
  410.     initial_user_parameter_temp[index_v] =
  411.     OPTIONS->INITIAL_PARAMETER_TEMPERATURE;
  412.     }
  413.  
  414.   if ((temperature_scale_parameters =
  415.        (double *) calloc (*number_parameters, sizeof (double))
  416.       ) == NULL)
  417.       exit (9);
  418.  
  419.   if (OPTIONS->USER_INITIAL_COST_TEMP == TRUE)
  420.     *initial_cost_temperature = *current_cost_temperature =
  421.       OPTIONS->user_cost_temperature[0];
  422.  
  423.   /* set parameters to the initial parameter values */
  424.   VFOR (index_v)
  425.     last_saved_state->parameter[index_v] =
  426.     current_generated_state->parameter[index_v] =
  427.     parameter_initial_final[index_v];
  428.  
  429.   /* test OPTIONS->SEQUENTIAL_PARAMETERS */
  430.   if (OPTIONS->SEQUENTIAL_PARAMETERS >= *number_parameters)
  431.     {
  432. #if ASA_PRINT
  433.       fprintf (ptr_asa_out, "OPTIONS->SEQUENTIAL_PARAMETERS > dimension-1\n");
  434.       fflush (ptr_asa_out);
  435. #endif
  436.       exit (9);
  437.     }
  438.  
  439.   /* save initial user value of OPTIONS->SEQUENTIAL_PARAMETERS */
  440.   *start_sequence = OPTIONS->SEQUENTIAL_PARAMETERS;
  441.  
  442. #if ASA_PRINT
  443. #if INT_ALLOC
  444.   fprintf (ptr_asa_out,
  445.        "*number_parameters = %d\n\n", *number_parameters);
  446. #else
  447.   fprintf (ptr_asa_out,
  448. #if INT_LONG
  449.        "*number_parameters = %ld\n\n", *number_parameters);
  450. #else
  451.        "*number_parameters = %d\n\n", *number_parameters);
  452. #endif
  453. #endif
  454.  
  455.   /* print the min, max, current values, and types of parameters */
  456.   fprintf (ptr_asa_out,
  457.        "index_v parameter_minimum parameter_maximum\
  458.  parameter_value parameter_type \n");
  459.  
  460. #if ASA_PRINT_INTERMED
  461.   VFOR (index_v)
  462. #if INT_ALLOC
  463.     fprintf (ptr_asa_out, " %-8d %-18.7g %-17.7g %-15.7g %-7d\n",
  464. #else
  465. #if INT_LONG
  466.     fprintf (ptr_asa_out, " %-8ld %-18.7g %-17.7g %-15.7g %-7d\n",
  467. #else
  468.     fprintf (ptr_asa_out, " %-8d %-18.7g %-17.7g %-15.7g %-7d\n",
  469. #endif
  470. #endif
  471.          index_v,
  472.          parameter_minimum[index_v],
  473.          parameter_maximum[index_v],
  474.          current_generated_state->parameter[index_v],
  475.          parameter_type[index_v]);
  476.  
  477.   fprintf (ptr_asa_out, "\n\n");
  478. #endif
  479.   /* Print out user-defined OPTIONS */
  480.  
  481.   if (OPTIONS->DELTA_PARAMETERS)
  482.     VFOR (index_v)
  483.       fprintf (ptr_asa_out,
  484. #if INT_ALLOC
  485.            "OPTIONS->user_delta_parameter[%d] = %12.7g\n",
  486. #else
  487. #if INT_LONG
  488.            "OPTIONS->user_delta_parameter[%ld] = %12.7g\n",
  489. #else
  490.            "OPTIONS->user_delta_parameter[%d] = %12.7g\n",
  491. #endif
  492. #endif
  493.            index_v, delta_parameter[index_v]);
  494.  
  495.   if (OPTIONS->QUENCH_PARAMETERS)
  496.     VFOR (index_v)
  497.       fprintf (ptr_asa_out,
  498. #if INT_ALLOC
  499.            "OPTIONS->user_quench_param_scale[%d] = %12.7g\n",
  500. #else
  501. #if INT_LONG
  502.            "OPTIONS->user_quench_param_scale[%ld] = %12.7g\n",
  503. #else
  504.            "OPTIONS->user_quench_param_scale[%d] = %12.7g\n",
  505. #endif
  506. #endif
  507.            index_v, quench_param[index_v]);
  508.  
  509.   if (OPTIONS->QUENCH_COST)
  510.     fprintf (ptr_asa_out,
  511.          "OPTIONS->user_quench_cost_scale = %12.7g\n",
  512.          *quench_cost);
  513.  
  514.   if (OPTIONS->USER_INITIAL_PARAMETERS_TEMPS)
  515.     VFOR (index_v)
  516.       fprintf (ptr_asa_out,
  517. #if INT_ALLOC
  518.            "OPTIONS->user_parameter_temperature[%d] = %12.7g\n",
  519. #else
  520. #if INT_LONG
  521.            "OPTIONS->user_parameter_temperature[%ld] = %12.7g\n",
  522. #else
  523.            "OPTIONS->user_parameter_temperature[%d] = %12.7g\n",
  524. #endif
  525. #endif
  526.            index_v, initial_user_parameter_temp[index_v]);
  527.  
  528.   if (OPTIONS->RATIO_TEMPERATURE_SCALES)
  529.     VFOR (index_v)
  530.       fprintf (ptr_asa_out,
  531. #if INT_ALLOC
  532.            "OPTIONS->user_temperature_ratio[%d] = %12.7g\n",
  533. #else
  534. #if INT_LONG
  535.            "OPTIONS->user_temperature_ratio[%ld] = %12.7g\n",
  536. #else
  537.            "OPTIONS->user_temperature_ratio[%d] = %12.7g\n",
  538. #endif
  539. #endif
  540.            index_v, OPTIONS->user_temperature_ratio[index_v]);
  541.  
  542.   if (OPTIONS->USER_INITIAL_COST_TEMP)
  543.     fprintf (ptr_asa_out,
  544.          "OPTIONS->user_cost_temperature[0] = %12.7g\n",
  545.          *initial_cost_temperature);
  546.  
  547.   fflush (ptr_asa_out);
  548. #endif
  549.  
  550. #if ASA_PARALLEL
  551. #if ASA_PRINT
  552.   fprintf (ptr_asa_out,
  553. #if INT_LONG
  554.        "Initial ASA_PARALLEL OPTIONS->\n\t gener_block = %ld\n\
  555.  \t gener_block_max = %ld\n \t gener_mov_avr= %d\n\n",
  556. #else
  557.        "ASA_PARALLEL OPTIONS->\n\t gener_block = %d\n\
  558.  \t gener_block_max = %d\n \t gener_mov_avr= %d\n\n",
  559. #endif
  560.     OPTIONS->gener_block, OPTIONS->gener_block_max, OPTIONS->gener_mov_avr);
  561. #endif
  562. #endif
  563.  
  564.   /* set the initial index of parameter generations to 1 */
  565.   VFOR (index_v)
  566.     index_parameter_generations[index_v] = 1;
  567.  
  568.   if (OPTIONS->USER_INITIAL_COST_TEMP == FALSE)
  569.     {
  570.       /* calculate the average cost over samplings of the cost function */
  571.       sampled_cost_sum = ZERO;
  572.       for (index_cost_constraint = 0;
  573.        index_cost_constraint < (OPTIONS->NUMBER_COST_SAMPLES);
  574.        ++index_cost_constraint)
  575.     {
  576.       *number_invalid_generated_states = 0;
  577.       repeated_invalid_states = 0;
  578.       OPTIONS->SEQUENTIAL_PARAMETERS = *start_sequence - 1;
  579.       do
  580.         {
  581.           ++(*number_invalid_generated_states);
  582.           ++repeated_invalid_states;
  583.           generate_new_state (user_random_generator,
  584.                   parameter_minimum,
  585.                   parameter_maximum,
  586.                   current_user_parameter_temp,
  587.                   number_parameters,
  588.                   parameter_type,
  589.                   current_generated_state,
  590.                   last_saved_state,
  591.                   OPTIONS);
  592.           *valid_state_generated_flag = TRUE;
  593.           current_generated_state->cost =
  594.         user_cost_function (current_generated_state->parameter,
  595.                     parameter_minimum,
  596.                     parameter_maximum,
  597.                     tangents,
  598.                     curvature,
  599.                     number_parameters,
  600.                     parameter_type,
  601.                     valid_state_generated_flag,
  602.                     exit_status,
  603.                     OPTIONS);
  604.           if (repeated_invalid_states >
  605.           OPTIONS->LIMIT_INVALID_GENERATED_STATES)
  606.         {
  607.           *exit_status = (int) TOO_MANY_INVALID_STATES;
  608.           goto EXIT_ASA;
  609.         }
  610.         }
  611.       while (*valid_state_generated_flag == FALSE);
  612.       --(*number_invalid_generated_states);
  613.       sampled_cost_sum += fabs (current_generated_state->cost);
  614.     }
  615.       /* set all cost temperatures to the average cost */
  616.       *initial_cost_temperature = *current_cost_temperature =
  617.     sampled_cost_sum / (OPTIONS->NUMBER_COST_SAMPLES);
  618.     }
  619.   /* set all parameters to the initial parameter values */
  620.   VFOR (index_v)
  621.     best_generated_state->parameter[index_v] =
  622.     last_saved_state->parameter[index_v] =
  623.     current_generated_state->parameter[index_v] =
  624.     parameter_initial_final[index_v];
  625.  
  626.   /* if using user's initial parameters */
  627.   if (OPTIONS->USER_INITIAL_PARAMETERS == TRUE)
  628.     {
  629.       *valid_state_generated_flag = TRUE;
  630.       current_generated_state->cost =
  631.     user_cost_function (current_generated_state->parameter,
  632.                 parameter_minimum,
  633.                 parameter_maximum,
  634.                 tangents,
  635.                 curvature,
  636.                 number_parameters,
  637.                 parameter_type,
  638.                 valid_state_generated_flag,
  639.                 exit_status,
  640.                 OPTIONS);
  641. #if ASA_PRINT
  642.       if (*valid_state_generated_flag == FALSE)
  643.     fprintf (ptr_asa_out,
  644.          "user's initial parameters generated \
  645. FALSE *valid_state_generated_flag\n");
  646. #endif
  647.     }
  648.   else
  649.     {
  650.       /* let asa generate valid initial parameters */
  651.       repeated_invalid_states = 0;
  652.       OPTIONS->SEQUENTIAL_PARAMETERS = *start_sequence - 1;
  653.       do
  654.     {
  655.       ++(*number_invalid_generated_states);
  656.       ++repeated_invalid_states;
  657.       generate_new_state (user_random_generator,
  658.                   parameter_minimum,
  659.                   parameter_maximum,
  660.                   current_user_parameter_temp,
  661.                   number_parameters,
  662.                   parameter_type,
  663.                   current_generated_state,
  664.                   last_saved_state,
  665.                   OPTIONS);
  666.       *valid_state_generated_flag = TRUE;
  667.       current_generated_state->cost =
  668.         user_cost_function (current_generated_state->parameter,
  669.                 parameter_minimum,
  670.                 parameter_maximum,
  671.                 tangents,
  672.                 curvature,
  673.                 number_parameters,
  674.                 parameter_type,
  675.                 valid_state_generated_flag,
  676.                 exit_status,
  677.                 OPTIONS);
  678.       if (repeated_invalid_states >
  679.           OPTIONS->LIMIT_INVALID_GENERATED_STATES)
  680.         {
  681.           *exit_status = (int) TOO_MANY_INVALID_STATES;
  682.           goto EXIT_ASA;
  683.         }
  684.     }
  685.       while (*valid_state_generated_flag == FALSE);
  686.       --(*number_invalid_generated_states);
  687.     }                /* OPTIONS->USER_INITIAL_PARAMETERS */
  688.  
  689.   /* set all states to the last one generated */
  690.   VFOR (index_v)
  691.   {
  692.     /* ignore parameters that have too small a range */
  693.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  694.       continue;
  695.     best_generated_state->parameter[index_v] =
  696.       last_saved_state->parameter[index_v] =
  697.       current_generated_state->parameter[index_v];
  698.   }
  699.  
  700.   /* set all costs to the last one generated */
  701.   best_generated_state->cost = last_saved_state->cost =
  702.     current_generated_state->cost;
  703.  
  704.   *accepted_to_generated_ratio = ONE;
  705.  
  706.   tmp_var_db1 = -log ((OPTIONS->TEMPERATURE_RATIO_SCALE));
  707.  
  708.   /* compute temperature scales */
  709.   tmp_var_db2 = log (OPTIONS->TEMPERATURE_ANNEAL_SCALE);
  710.   temperature_scale = tmp_var_db1 * exp (-tmp_var_db2
  711.                      / (double) *number_parameters);
  712.  
  713.   if (OPTIONS->RATIO_TEMPERATURE_SCALES)
  714.     {
  715.       VFOR (index_v)
  716.     temperature_scale_parameters[index_v] =
  717.     tmp_var_db1 * exp (-(tmp_var_db2 * quench_param[index_v])
  718.                / (double) *number_parameters)
  719.     * OPTIONS->user_temperature_ratio[index_v];
  720.     }
  721.   else
  722.     {
  723.       VFOR (index_v)
  724.     temperature_scale_parameters[index_v] =
  725.     tmp_var_db1 * exp (-(tmp_var_db2 * quench_param[index_v])
  726.                / (double) *number_parameters);
  727.     }
  728.  
  729.   *temperature_scale_cost =
  730.     tmp_var_db1 * exp (-(tmp_var_db2 * quench_cost[0])
  731.                / (double) *number_parameters)
  732.     * OPTIONS->COST_PARAMETER_SCALE;
  733.  
  734.   /* do not calculate curvatures initially */
  735.   *curvature_flag = FALSE;
  736.  
  737. #if ASA_PRINT
  738.   fprintf (ptr_asa_out,
  739.        "temperature_scale = %12.7g\n",
  740.        temperature_scale);
  741.   if (OPTIONS->RATIO_TEMPERATURE_SCALES)
  742.     {
  743. #if ASA_PRINT_INTERMED
  744.       VFOR (index_v)
  745.       {
  746.     fprintf (ptr_asa_out,
  747. #if INT_ALLOC
  748.          "temperature_scale_parameters[%d] = %12.7g\n",
  749. #else
  750. #if INT_LONG
  751.          "temperature_scale_parameters[%ld] = %12.7g\n",
  752. #else
  753.          "temperature_scale_parameters[%d] = %12.7g\n",
  754. #endif
  755. #endif
  756.          index_v, temperature_scale_parameters[index_v]);
  757.       }
  758. #endif
  759.     }
  760.   else
  761.     {
  762.       fprintf (ptr_asa_out,
  763.            "temperature_scale_parameters[0] = %12.7g\n",
  764.            temperature_scale_parameters[0]);
  765.     }
  766.   fprintf (ptr_asa_out,
  767.        "*temperature_scale_cost = %12.7g\n",
  768.        *temperature_scale_cost);
  769.   fprintf (ptr_asa_out, "\n\n");
  770.  
  771. #if ASA_PRINT_INTERMED
  772.   print_state (parameter_minimum,
  773.            parameter_maximum,
  774.            tangents,
  775.            curvature,
  776.            current_cost_temperature,
  777.            current_user_parameter_temp,
  778.            accepted_to_generated_ratio,
  779.            number_parameters,
  780.            curvature_flag,
  781.            number_accepted,
  782.            index_cost_acceptances,
  783.            number_generated,
  784.            number_invalid_generated_states,
  785.            last_saved_state,
  786.            best_generated_state,
  787.            ptr_asa_out,
  788.            OPTIONS);
  789. #endif
  790.   fprintf (ptr_asa_out, "\n");
  791.  
  792.   fflush (ptr_asa_out);
  793. #endif
  794.  
  795. #if ASA_SAMPLE
  796. #if ASA_PRINT
  797.   fprintf (ptr_asa_out,
  798.        ":SAMPLE:   n_accept   cost        cost_temp    bias_accept    \
  799.  aver_weight\n");
  800.   fprintf (ptr_asa_out,
  801.        ":SAMPLE:   index      param[]     temp[]       bias_gener[]   \
  802.  range[]\n");
  803. #endif
  804. #endif
  805.  
  806.   /* reset the current cost and the number of generations performed */
  807.   *number_invalid_generated_states = 0;
  808.   *best_number_generated_saved =
  809.     *number_generated = *recent_number_generated = 0;
  810.   VFOR (index_v)
  811.   {
  812.     /* ignore parameters that have too small a range */
  813.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  814.       continue;
  815.     index_parameter_generations[index_v] = 1;
  816.   }
  817.  
  818.   OPTIONS->SEQUENTIAL_PARAMETERS = *start_sequence - 1;
  819.  
  820.   /* MAIN ANNEALING LOOP */
  821.   while (
  822.       ((*number_accepted <= OPTIONS->LIMIT_ACCEPTANCES)
  823.        || (OPTIONS->LIMIT_ACCEPTANCES == 0))
  824.       &&
  825.       ((*number_generated <= OPTIONS->LIMIT_GENERATED)
  826.        || (OPTIONS->LIMIT_GENERATED == 0))
  827.     )
  828.     {
  829.  
  830.       tmp_var_db1 = -log ((OPTIONS->TEMPERATURE_RATIO_SCALE));
  831.  
  832.       /* compute temperature scales */
  833.       tmp_var_db2 = log (OPTIONS->TEMPERATURE_ANNEAL_SCALE);
  834.       temperature_scale = tmp_var_db1 * exp (-tmp_var_db2
  835.                          / (double) *number_parameters);
  836.  
  837.       if (OPTIONS->RATIO_TEMPERATURE_SCALES)
  838.     {
  839.       VFOR (index_v)
  840.         temperature_scale_parameters[index_v] =
  841.         tmp_var_db1 * exp (-(tmp_var_db2 * quench_param[index_v])
  842.                    / (double) *number_parameters)
  843.         * OPTIONS->user_temperature_ratio[index_v];
  844.     }
  845.       else
  846.     {
  847.       VFOR (index_v)
  848.         temperature_scale_parameters[index_v] =
  849.         tmp_var_db1 * exp (-(tmp_var_db2 * quench_param[index_v])
  850.                    / (double) *number_parameters);
  851.     }
  852.  
  853.       *temperature_scale_cost =
  854.     tmp_var_db1 * exp (-(tmp_var_db2 * quench_cost[0])
  855.                / (double) *number_parameters)
  856.     * OPTIONS->COST_PARAMETER_SCALE;
  857.  
  858.       /* CALCULATE NEW TEMPERATURES */
  859.  
  860.       /* calculate new parameter temperatures */
  861.       VFOR (index_v)
  862.       {
  863.     /* skip parameters with too small a range */
  864.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  865.       continue;
  866.  
  867.     log_new_temperature_ratio =
  868.       -temperature_scale_parameters[index_v] *
  869.       pow ((double) index_parameter_generations[index_v],
  870.            quench_param[index_v] / (double) *number_parameters);
  871.     /* check (and correct) for too large an exponent */
  872.     log_new_temperature_ratio =
  873.       EXPONENT_CHECK (log_new_temperature_ratio);
  874.     current_user_parameter_temp[index_v] =
  875.       initial_user_parameter_temp[index_v]
  876.       * exp (log_new_temperature_ratio);
  877.  
  878. #if NO_PARAM_TEMP_TEST
  879.     if (current_user_parameter_temp[index_v]
  880.         < (double) EPS_DOUBLE)
  881.       current_user_parameter_temp[index_v] =
  882.         (double) EPS_DOUBLE;
  883. #else
  884.     /* check for too small a parameter temperature */
  885.     if (current_user_parameter_temp[index_v]
  886.         < (double) EPS_DOUBLE)
  887.       {
  888.         *exit_status = (int) P_TEMP_TOO_SMALL;
  889.         *index_exit_v = index_v;
  890.         goto EXIT_ASA;
  891.       }
  892. #endif
  893.       }
  894.  
  895.       /* calculate new cost temperature */
  896.       log_new_temperature_ratio =
  897.     -*temperature_scale_cost
  898.     * pow ((double) *index_cost_acceptances,
  899.            (*quench_cost) / (double) *number_parameters);
  900.       log_new_temperature_ratio =
  901.     EXPONENT_CHECK (log_new_temperature_ratio);
  902.       *current_cost_temperature = *initial_cost_temperature
  903.     * exp (log_new_temperature_ratio);
  904.  
  905. #if NO_COST_TEMP_TEST
  906.       if (*current_cost_temperature < (double) EPS_DOUBLE)
  907.     *current_cost_temperature = (double) EPS_DOUBLE;
  908. #else
  909.       /* check for too small a cost temperature */
  910.       if (*current_cost_temperature < (double) EPS_DOUBLE)
  911.     {
  912.       *exit_status = (int) C_TEMP_TOO_SMALL;
  913.       goto EXIT_ASA;
  914.     }
  915. #endif
  916.       /* GENERATE NEW PARAMETERS */
  917.  
  918.       /* generate a new valid set of parameters */
  919. #if ASA_PARALLEL
  920.       /* *** ENTER CODE TO SPAWN OFF PARALLEL GENERATED STATES *** */
  921.       if (OPTIONS->gener_block_max > parallel_block_max)
  922.     {
  923.       for (index_parallel = 0; index_parallel < parallel_block_max;
  924.            ++index_parallel)
  925.         {
  926.           free (gener_block_state[index_parallel].parameter);
  927.         }
  928.       free (gener_block_state);
  929.  
  930.       if ((gener_block_state =
  931.            (STATE *) calloc (OPTIONS->gener_block_max,
  932.                  sizeof (STATE))) == NULL)
  933.         exit (9);
  934.  
  935.       parallel_block_max = OPTIONS->gener_block_max;
  936.  
  937.       for (index_parallel = 0; index_parallel < parallel_block_max;
  938.            ++index_parallel)
  939.         {
  940.           if ((gener_block_state[index_parallel].parameter =
  941.            (double *) calloc (*number_parameters, sizeof (double))
  942.           ) == NULL)
  943.           exit (9);
  944.         }
  945.     }
  946.       for (index_parallel = 0; index_parallel < OPTIONS->gener_block;
  947.        ++index_parallel)
  948.     {
  949. #endif /* ASA_PARALLEL */
  950.  
  951.       repeated_invalid_states = 0;
  952.       do
  953.         {
  954.           ++(*number_invalid_generated_states);
  955.           ++repeated_invalid_states;
  956.           generate_new_state (user_random_generator,
  957.                   parameter_minimum,
  958.                   parameter_maximum,
  959.                   current_user_parameter_temp,
  960.                   number_parameters,
  961.                   parameter_type,
  962.                   current_generated_state,
  963.                   last_saved_state,
  964.                   OPTIONS);
  965.           *valid_state_generated_flag = TRUE;
  966.           current_generated_state->cost =
  967.         user_cost_function (current_generated_state->parameter,
  968.                     parameter_minimum,
  969.                     parameter_maximum,
  970.                     tangents,
  971.                     curvature,
  972.                     number_parameters,
  973.                     parameter_type,
  974.                     valid_state_generated_flag,
  975.                     exit_status,
  976.                     OPTIONS);
  977.           if (repeated_invalid_states >
  978.           OPTIONS->LIMIT_INVALID_GENERATED_STATES)
  979.         {
  980.           *exit_status = (int) TOO_MANY_INVALID_STATES;
  981.           goto EXIT_ASA;
  982.         }
  983.         }
  984.       while (*valid_state_generated_flag == FALSE);
  985.       --(*number_invalid_generated_states);
  986. #if ASA_PARALLEL
  987.       gener_block_state[index_parallel].cost
  988.         = current_generated_state->cost;
  989.       VFOR (index_v)
  990.       {
  991.         /* ignore parameters with too small a range */
  992.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  993.           continue;
  994.         gener_block_state[index_parallel].parameter[index_v] =
  995.           current_generated_state->parameter[index_v];
  996.       }
  997.     }
  998.       /* *** EXIT CODE SPAWNING OFF PARALLEL GENERATED STATES *** */
  999. #endif /* ASA_PARALLEL */
  1000.  
  1001.       /* ACCEPT/REJECT NEW PARAMETERS */
  1002.  
  1003. #if ASA_PARALLEL
  1004.       for (index_parallel = 0; index_parallel < OPTIONS->gener_block;
  1005.        ++index_parallel)
  1006.     {
  1007.       current_generated_state->cost
  1008.         = gener_block_state[index_parallel].cost;
  1009.       VFOR (index_v)
  1010.       {
  1011.         /* ignore parameters with too small a range */
  1012.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1013.           continue;
  1014.         current_generated_state->parameter[index_v] =
  1015.           gener_block_state[index_parallel].parameter[index_v];
  1016.       }
  1017. #endif /* ASA_PARALLEL */
  1018.       /* decide to accept/reject the new state */
  1019.       accept_new_state (user_random_generator,
  1020.                 parameter_minimum,
  1021.                 parameter_maximum,
  1022.                 current_cost_temperature,
  1023. #if ASA_SAMPLE
  1024.                 current_user_parameter_temp,
  1025. #endif
  1026.                 number_parameters,
  1027.                 recent_number_acceptances,
  1028.                 number_accepted,
  1029.                 index_cost_acceptances,
  1030.                 number_acceptances_saved,
  1031.                 recent_number_generated,
  1032.                 number_generated,
  1033.                 index_parameter_generations,
  1034.                 current_generated_state,
  1035.                 last_saved_state,
  1036. #if ASA_SAMPLE
  1037.                 ptr_asa_out,
  1038. #endif
  1039.                 OPTIONS);
  1040.  
  1041.       /* calculate the ratio of acceptances to generated states */
  1042.       *accepted_to_generated_ratio =
  1043.         (double) (*recent_number_acceptances + 1) /
  1044.         (double) (*recent_number_generated + 1);
  1045.  
  1046.       /* CHECK FOR NEW MINIMUM */
  1047.  
  1048.       if (current_generated_state->cost < best_generated_state->cost)
  1049.         {
  1050.           /* NEW MINIMUM FOUND */
  1051.  
  1052.           /* reset the recent acceptances and generated counts */
  1053. #if ASA_PARALLEL
  1054.           parallel_generated = *recent_number_generated;
  1055. #endif
  1056.           *recent_number_acceptances = *recent_number_generated = 0;
  1057.           *best_number_generated_saved = *number_generated;
  1058.           *best_number_accepted_saved = *number_accepted;
  1059.           index_cost_repeat = 0;
  1060.  
  1061.           /* copy the current state into the best_generated state */
  1062.           best_generated_state->cost = current_generated_state->cost;
  1063.           VFOR (index_v)
  1064.           {
  1065.         /* ignore parameters that have too small a range */
  1066.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1067.           continue;
  1068.         best_generated_state->parameter[index_v] =
  1069.           current_generated_state->parameter[index_v];
  1070.           }
  1071.  
  1072.           /* printout the new minimum state and value */
  1073. #if ASA_PRINT
  1074. #if INT_LONG
  1075.           fprintf (ptr_asa_out,
  1076.                "best...->cost=%-12.7g  \
  1077. *number_accepted=%ld  *number_generated=%ld\n",
  1078.                best_generated_state->cost,
  1079.                *number_accepted, *number_generated);
  1080. #else
  1081.           fprintf (ptr_asa_out,
  1082.                "best...->cost=%-12.7g  \
  1083. *number_accepted=%d  *number_generated=%d\n",
  1084.                best_generated_state->cost,
  1085.                *number_accepted, *number_generated);
  1086. #endif /* INT_LOG */
  1087. #if ASA_PARALLEL
  1088.           /* print OPTIONS->gener_block just used */
  1089.           fprintf (ptr_asa_out,
  1090. #if INT_LONG
  1091.                "OPTIONS->gener_block = %ld\n", OPTIONS->gener_block);
  1092. #else
  1093.                "OPTIONS->gener_block = %d\n", OPTIONS->gener_block);
  1094. #endif /* INT_LONG */
  1095. #endif /* ASA_PARALLEL */
  1096. #if ASA_PRINT_MORE
  1097.           print_state (parameter_minimum,
  1098.                parameter_maximum,
  1099.                tangents,
  1100.                curvature,
  1101.                current_cost_temperature,
  1102.                current_user_parameter_temp,
  1103.                accepted_to_generated_ratio,
  1104.                number_parameters,
  1105.                curvature_flag,
  1106.                number_accepted,
  1107.                index_cost_acceptances,
  1108.                number_generated,
  1109.                number_invalid_generated_states,
  1110.                last_saved_state,
  1111.                best_generated_state,
  1112.                ptr_asa_out,
  1113.                OPTIONS);
  1114. #endif /* ASA_PRINT_MORE */
  1115.           fflush (ptr_asa_out);
  1116. #endif /* ASA_PRINT */
  1117. #if ASA_PARALLEL
  1118.           /* leave index_parallel loop after new minimum */
  1119.           break;
  1120. #endif /* ASA_PARALLEL */
  1121.         }
  1122. #if ASA_PARALLEL
  1123.     }
  1124. #endif /* ASA_PARALLEL */
  1125. #if ASA_PARALLEL
  1126.       if (OPTIONS->gener_mov_avr > 0)
  1127.     {
  1128.       OPTIONS->gener_block = (LONG_INT)
  1129.         ((((double) OPTIONS->gener_mov_avr - ONE)
  1130.         * (double) (OPTIONS->gener_block) + (double) parallel_generated)
  1131.          / (double) (OPTIONS->gener_mov_avr));
  1132.       OPTIONS->gener_block =
  1133.         MIN (OPTIONS->gener_block, parallel_block_max);
  1134.     }
  1135. #endif /* ASA_PARALLEL */
  1136.  
  1137.       /* PERIODIC TESTING/REANNEALING/PRINTING SECTION */
  1138.  
  1139.       tmp_var_int = (int) (*number_accepted %
  1140.                ((LONG_INT) OPTIONS->TESTING_FREQUENCY_MODULUS));
  1141.       if ((tmp_var_int == 0 && *number_acceptances_saved
  1142.        == *number_accepted)
  1143.       || *accepted_to_generated_ratio
  1144.       < (OPTIONS->ACCEPTED_TO_GENERATED_RATIO))
  1145.     {
  1146.       if (*accepted_to_generated_ratio
  1147.           < (OPTIONS->ACCEPTED_TO_GENERATED_RATIO))
  1148.         *recent_number_acceptances = *recent_number_generated = 0;
  1149.  
  1150.       /* if best.cost repeats OPTIONS->MAXIMUM_COST_REPEAT then exit */
  1151.       if (fabs (last_saved_state->cost
  1152.             - best_generated_state->cost) <
  1153.           OPTIONS->COST_PRECISION)
  1154.         {
  1155.           ++index_cost_repeat;
  1156.           if (index_cost_repeat == (OPTIONS->MAXIMUM_COST_REPEAT))
  1157.         {
  1158.           *exit_status = (int) COST_REPEATING;
  1159.           goto EXIT_ASA;
  1160.         }
  1161.         }
  1162.       else
  1163.         {
  1164.           index_cost_repeat = 0;
  1165.         }
  1166.  
  1167.       /* set to FALSE to skip reannealing */
  1168.       if (OPTIONS->ACTIVATE_REANNEAL == TRUE)
  1169.         {
  1170.           /* calculate tangents, not curvatures, to reanneal */
  1171.           *curvature_flag = FALSE;
  1172.           cost_derivatives (user_cost_function,
  1173.                 parameter_minimum,
  1174.                 parameter_maximum,
  1175.                 tangents,
  1176.                 curvature,
  1177.                 maximum_tangent,
  1178.                 delta_parameter,
  1179.                 number_parameters,
  1180.                 parameter_type,
  1181.                 exit_status,
  1182.                 curvature_flag,
  1183.                 valid_state_generated_flag,
  1184.                 number_invalid_generated_states,
  1185.                 current_generated_state,
  1186.                 best_generated_state,
  1187.                 ptr_asa_out,
  1188.                 OPTIONS);
  1189.           reanneal (parameter_minimum,
  1190.             parameter_maximum,
  1191.             tangents,
  1192.             maximum_tangent,
  1193.             current_cost_temperature,
  1194.             initial_cost_temperature,
  1195.             temperature_scale_cost,
  1196.             current_user_parameter_temp,
  1197.             initial_user_parameter_temp,
  1198.             temperature_scale_parameters,
  1199.             quench_param,
  1200.             quench_cost,
  1201.             number_parameters,
  1202.             parameter_type,
  1203.             index_cost_acceptances,
  1204.             index_parameter_generations,
  1205.             last_saved_state,
  1206.             best_generated_state,
  1207.             OPTIONS);
  1208.         }
  1209. #if ASA_PRINT_INTERMED
  1210. #if ASA_PRINT
  1211.       print_state (parameter_minimum,
  1212.                parameter_maximum,
  1213.                tangents,
  1214.                curvature,
  1215.                current_cost_temperature,
  1216.                current_user_parameter_temp,
  1217.                accepted_to_generated_ratio,
  1218.                number_parameters,
  1219.                curvature_flag,
  1220.                number_accepted,
  1221.                index_cost_acceptances,
  1222.                number_generated,
  1223.                number_invalid_generated_states,
  1224.                last_saved_state,
  1225.                best_generated_state,
  1226.                ptr_asa_out,
  1227.                OPTIONS);
  1228.  
  1229.       fprintf (ptr_asa_out, "\n");
  1230.       fflush (ptr_asa_out);
  1231. #endif
  1232. #endif
  1233.     }
  1234.     }
  1235.  
  1236.   /* FINISHED ANNEALING and MINIMIZATION */
  1237.  
  1238.   *exit_status = (int) NORMAL_EXIT;
  1239. EXIT_ASA:
  1240.   asa_exit (user_cost_function,
  1241.         &final_cost,
  1242.         parameter_initial_final,
  1243.         parameter_minimum,
  1244.         parameter_maximum,
  1245.         tangents,
  1246.         curvature,
  1247.         maximum_tangent,
  1248.         delta_parameter,
  1249.         current_cost_temperature,
  1250.         initial_user_parameter_temp,
  1251.         current_user_parameter_temp,
  1252.         accepted_to_generated_ratio,
  1253.         number_parameters,
  1254.         parameter_type,
  1255.         valid_state_generated_flag,
  1256.         exit_status,
  1257.         index_exit_v,
  1258.         start_sequence,
  1259.         number_accepted,
  1260.         best_number_accepted_saved,
  1261.         index_cost_acceptances,
  1262.         number_generated,
  1263.         number_invalid_generated_states,
  1264.         index_parameter_generations,
  1265.         best_number_generated_saved,
  1266.         current_generated_state,
  1267.         last_saved_state,
  1268.         best_generated_state,
  1269.         ptr_asa_out,
  1270.         OPTIONS);
  1271.   free (curvature_flag);
  1272.   free (maximum_tangent);
  1273.   free (accepted_to_generated_ratio);
  1274.   free (temperature_scale_cost);
  1275.   free (current_cost_temperature);
  1276.   free (initial_cost_temperature);
  1277.   free (number_generated);
  1278.   free (best_number_generated_saved);
  1279.   free (recent_number_generated);
  1280.   free (number_accepted);
  1281.   free (recent_number_acceptances);
  1282.   free (index_cost_acceptances);
  1283.   free (number_acceptances_saved);
  1284.   free (best_number_accepted_saved);
  1285.   free (number_invalid_generated_states);
  1286.   free (current_generated_state->parameter);
  1287.   free (last_saved_state->parameter);
  1288.   free (best_generated_state->parameter);
  1289.   free (current_generated_state);
  1290.   free (last_saved_state);
  1291.   free (best_generated_state);
  1292. #if ASA_PARALLEL
  1293.   for (index_parallel = 0; index_parallel < parallel_block_max;
  1294.        ++index_parallel)
  1295.     {
  1296.       free (gener_block_state[index_parallel].parameter);
  1297.     }
  1298.   free (gener_block_state);
  1299. #endif
  1300.   free (initial_user_parameter_temp);
  1301.   free (index_exit_v);
  1302.   free (start_sequence);
  1303.   free (index_parameter_generations);
  1304.   if (OPTIONS->DELTA_PARAMETERS == FALSE)
  1305.     free (delta_parameter);
  1306.   if (OPTIONS->QUENCH_PARAMETERS == FALSE)
  1307.     free (quench_param);
  1308.   if (OPTIONS->QUENCH_COST == FALSE)
  1309.     free (quench_cost);
  1310.   free (current_user_parameter_temp);
  1311.   free (temperature_scale_parameters);
  1312.   if (recursive_asa_open == 0)
  1313.     asa_open = FALSE;
  1314.   return (final_cost);
  1315. }
  1316.  
  1317. /***********************************************************************
  1318. * asa_exit
  1319. *    This procedures copies the best parameters and cost into
  1320. *       final_cost and parameter_initial_final
  1321. ***********************************************************************/
  1322. #if HAVE_ANSI
  1323. void
  1324. asa_exit (double (*user_cost_function) (
  1325.                double *, double *, double *, double *, double *,
  1326.               ALLOC_INT *, int *, int *, int *, USER_DEFINES *),
  1327.       double *final_cost,
  1328.       double *parameter_initial_final,
  1329.       double *parameter_minimum,
  1330.       double *parameter_maximum,
  1331.       double *tangents,
  1332.       double *curvature,
  1333.       double *maximum_tangent,
  1334.       double *delta_parameter,
  1335.       double *current_cost_temperature,
  1336.       double *initial_user_parameter_temp,
  1337.       double *current_user_parameter_temp,
  1338.       double *accepted_to_generated_ratio,
  1339.       ALLOC_INT * number_parameters,
  1340.       int *parameter_type,
  1341.       int *valid_state_generated_flag,
  1342.       int *exit_status,
  1343.       ALLOC_INT * index_exit_v,
  1344.       ALLOC_INT * start_sequence,
  1345.       LONG_INT * number_accepted,
  1346.       LONG_INT * best_number_accepted_saved,
  1347.       LONG_INT * index_cost_acceptances,
  1348.       LONG_INT * number_generated,
  1349.       LONG_INT * number_invalid_generated_states,
  1350.       LONG_INT * index_parameter_generations,
  1351.       LONG_INT * best_number_generated_saved,
  1352.       STATE * current_generated_state,
  1353.       STATE * last_saved_state,
  1354.       STATE * best_generated_state,
  1355.       FILE * ptr_asa_out,
  1356.       USER_DEFINES * OPTIONS)
  1357. #else
  1358. void
  1359. asa_exit (user_cost_function,
  1360.       final_cost,
  1361.       parameter_initial_final,
  1362.       parameter_minimum,
  1363.       parameter_maximum,
  1364.       tangents,
  1365.       curvature,
  1366.       maximum_tangent,
  1367.       delta_parameter,
  1368.       current_cost_temperature,
  1369.       initial_user_parameter_temp,
  1370.       current_user_parameter_temp,
  1371.       accepted_to_generated_ratio,
  1372.       number_parameters,
  1373.       parameter_type,
  1374.       valid_state_generated_flag,
  1375.       exit_status,
  1376.       index_exit_v,
  1377.       start_sequence,
  1378.       number_accepted,
  1379.       best_number_accepted_saved,
  1380.       index_cost_acceptances,
  1381.       number_generated,
  1382.       number_invalid_generated_states,
  1383.       index_parameter_generations,
  1384.       best_number_generated_saved,
  1385.       current_generated_state,
  1386.       last_saved_state,
  1387.       best_generated_state,
  1388.       ptr_asa_out,
  1389.       OPTIONS)
  1390.      double (*user_cost_function) ();
  1391.      double *final_cost;
  1392.      double *parameter_initial_final;
  1393.      double *parameter_minimum;
  1394.      double *parameter_maximum;
  1395.      double *tangents;
  1396.      double *curvature;
  1397.      double *maximum_tangent;
  1398.      double *delta_parameter;
  1399.      double *current_cost_temperature;
  1400.      double *initial_user_parameter_temp;
  1401.      double *current_user_parameter_temp;
  1402.      double *accepted_to_generated_ratio;
  1403.      ALLOC_INT *number_parameters;
  1404.      int *parameter_type;
  1405.      int *valid_state_generated_flag;
  1406.      int *exit_status;
  1407.      ALLOC_INT *index_exit_v;
  1408.      ALLOC_INT *start_sequence;
  1409.      LONG_INT *number_accepted;
  1410.      LONG_INT *best_number_accepted_saved;
  1411.      LONG_INT *index_cost_acceptances;
  1412.      LONG_INT *number_generated;
  1413.      LONG_INT *number_invalid_generated_states;
  1414.      LONG_INT *index_parameter_generations;
  1415.      LONG_INT *best_number_generated_saved;
  1416.      STATE *current_generated_state;
  1417.      STATE *last_saved_state;
  1418.      STATE *best_generated_state;
  1419.      FILE *ptr_asa_out;
  1420.      USER_DEFINES *OPTIONS;
  1421. #endif
  1422. {
  1423.   ALLOC_INT index_v;        /* iteration index */
  1424.   int *curvature_flag;
  1425.  
  1426.   if ((curvature_flag =
  1427.        (int *) calloc (1, sizeof (int))) == NULL)
  1428.       exit (9);
  1429.  
  1430.   if (*exit_status != TOO_MANY_INVALID_STATES)
  1431.     {
  1432.       /* calculate curvatures and tangents at best point */
  1433.       *curvature_flag = TRUE;
  1434.       cost_derivatives (user_cost_function,
  1435.             parameter_minimum,
  1436.             parameter_maximum,
  1437.             tangents,
  1438.             curvature,
  1439.             maximum_tangent,
  1440.             delta_parameter,
  1441.             number_parameters,
  1442.             parameter_type,
  1443.             exit_status,
  1444.             curvature_flag,
  1445.             valid_state_generated_flag,
  1446.             number_invalid_generated_states,
  1447.             current_generated_state,
  1448.             best_generated_state,
  1449.             ptr_asa_out,
  1450.             OPTIONS);
  1451.     }
  1452.   /* return final function minimum and associated parameters */
  1453.   *final_cost = best_generated_state->cost;
  1454.   VFOR (index_v)
  1455.   {
  1456.     parameter_initial_final[index_v] =
  1457.       best_generated_state->parameter[index_v];
  1458.   }
  1459.  
  1460. #if ASA_PRINT
  1461.   print_state (parameter_minimum,
  1462.            parameter_maximum,
  1463.            tangents,
  1464.            curvature,
  1465.            current_cost_temperature,
  1466.            current_user_parameter_temp,
  1467.            accepted_to_generated_ratio,
  1468.            number_parameters,
  1469.            curvature_flag,
  1470.            number_accepted,
  1471.            index_cost_acceptances,
  1472.            number_generated,
  1473.            number_invalid_generated_states,
  1474.            last_saved_state,
  1475.            best_generated_state,
  1476.            ptr_asa_out,
  1477.            OPTIONS);
  1478.  
  1479.   switch (*exit_status)
  1480.     {
  1481.     case NORMAL_EXIT:
  1482.       fprintf (ptr_asa_out,
  1483.            "\n\n NORMAL_EXIT exit_status = %d\n",
  1484.            *exit_status);
  1485.       break;
  1486.     case P_TEMP_TOO_SMALL:
  1487.       fprintf (ptr_asa_out,
  1488.            "\n\n P_TEMP_TOO_SMALL exit_status = %d\n",
  1489.            *exit_status);
  1490.       fprintf (ptr_asa_out,
  1491. #if INT_ALLOC
  1492.            "current_user_parameter_temp[%d] too small = %12.7g\n",
  1493. #else
  1494. #if INT_LONG
  1495.            "current_user_parameter_temp[%ld] too small = %12.7g\n",
  1496. #else
  1497.            "current_user_parameter_temp[%d] too small = %12.7g\n",
  1498. #endif
  1499. #endif
  1500.            *index_exit_v, current_user_parameter_temp[*index_exit_v]);
  1501.       break;
  1502.     case C_TEMP_TOO_SMALL:
  1503.       fprintf (ptr_asa_out,
  1504.            "\n\n C_TEMP_TOO_SMALL exit_status = %d\n",
  1505.            *exit_status);
  1506.       fprintf (ptr_asa_out,
  1507.            "*current_cost_temperature too small = %12.7g\n",
  1508.            *current_cost_temperature);
  1509.       break;
  1510.     case COST_REPEATING:
  1511.       fprintf (ptr_asa_out,
  1512.            "\n\n COST_REPEATING exit_status = %d\n",
  1513.            *exit_status);
  1514.       break;
  1515.     case TOO_MANY_INVALID_STATES:
  1516.       fprintf (ptr_asa_out,
  1517.            "\n\n  TOO_MANY_INVALID_STATES exit_status = %d\n",
  1518.            *exit_status);
  1519.       break;
  1520.     default:
  1521.       fprintf (ptr_asa_out, "\n\n ERR: no exit code available = %d\n",
  1522.            *exit_status);
  1523.     }
  1524.  
  1525.   fprintf (ptr_asa_out,
  1526.        "final_cost = best_generated_state->cost = %-12.7g\n",
  1527.        *final_cost);
  1528. #if INT_LONG
  1529.   fprintf (ptr_asa_out,
  1530.        "*number_accepted at best_generated_state->cost = %ld\n",
  1531.        *best_number_accepted_saved);
  1532.   fprintf (ptr_asa_out,
  1533.        "*number_generated at best_generated_state->cost = %ld\n",
  1534.        *best_number_generated_saved);
  1535. #else
  1536.   fprintf (ptr_asa_out,
  1537.        "*number_accepted at best_generated_state->cost = %d\n",
  1538.        *best_number_accepted_saved);
  1539.   fprintf (ptr_asa_out,
  1540.        "*number_generated at best_generated_state->cost = %d\n",
  1541.        *best_number_generated_saved);
  1542. #endif
  1543. #endif
  1544.  
  1545. #if ASA_TEMPLATE
  1546. #if OPTIONAL_DATA
  1547.   if (OPTIONS->asa_data[0] > (double) MIN_DOUBLE)
  1548.     OPTIONS->asa_data[1] = (double) (*best_number_generated_saved);
  1549. #endif
  1550. #endif
  1551.  
  1552.   /* reset OPTIONS->SEQUENTIAL_PARAMETERS */
  1553.   OPTIONS->SEQUENTIAL_PARAMETERS = *start_sequence;
  1554.  
  1555.   /* return unused space */
  1556.   free (curvature_flag);
  1557.  
  1558. #if ASA_PRINT
  1559. #if TIME_CALC
  1560.   /* print ending time */
  1561.   print_time ("asa_end", ptr_asa_out);
  1562. #endif
  1563.   fprintf (ptr_asa_out, "\n\n\n");
  1564.   fflush (ptr_asa_out);
  1565.   fclose (ptr_asa_out);
  1566. #endif
  1567. }
  1568.  
  1569. /***********************************************************************
  1570. * generate_new_state
  1571. *       Generates a valid new state from the old state
  1572. ***********************************************************************/
  1573. #if HAVE_ANSI
  1574. void
  1575. generate_new_state (double (*user_random_generator) (),
  1576.             double *parameter_minimum,
  1577.             double *parameter_maximum,
  1578.             double *current_user_parameter_temp,
  1579.             ALLOC_INT * number_parameters,
  1580.             int *parameter_type,
  1581.             STATE * current_generated_state,
  1582.             STATE * last_saved_state,
  1583.             USER_DEFINES * OPTIONS)
  1584. #else
  1585. void
  1586. generate_new_state (user_random_generator,
  1587.             parameter_minimum,
  1588.             parameter_maximum,
  1589.             current_user_parameter_temp,
  1590.             number_parameters,
  1591.             parameter_type,
  1592.             current_generated_state,
  1593.             last_saved_state,
  1594.             OPTIONS)
  1595.      double (*user_random_generator) ();
  1596.      double *parameter_minimum;
  1597.      double *parameter_maximum;
  1598.      double *current_user_parameter_temp;
  1599.      ALLOC_INT *number_parameters;
  1600.      int *parameter_type;
  1601.      STATE *current_generated_state;
  1602.      STATE *last_saved_state;
  1603.      USER_DEFINES *OPTIONS;
  1604. #endif
  1605. {
  1606.   ALLOC_INT index_v;
  1607.   double x;
  1608.   double parameter_v, min_parameter_v, max_parameter_v, temperature_v,
  1609.     parameter_range_v;
  1610.  
  1611.   /* generate a new value for each parameter */
  1612.   VFOR (index_v)
  1613.   {
  1614.     if (OPTIONS->SEQUENTIAL_PARAMETERS >= -1)
  1615.       {
  1616.     ++OPTIONS->SEQUENTIAL_PARAMETERS;
  1617.     if (OPTIONS->SEQUENTIAL_PARAMETERS == *number_parameters)
  1618.       OPTIONS->SEQUENTIAL_PARAMETERS = 0;
  1619.     index_v = OPTIONS->SEQUENTIAL_PARAMETERS;
  1620.       }
  1621.     min_parameter_v = parameter_minimum[index_v];
  1622.     max_parameter_v = parameter_maximum[index_v];
  1623.     parameter_range_v = max_parameter_v - min_parameter_v;
  1624.  
  1625.     /* ignore parameters that have too small a range */
  1626.     if (fabs (parameter_range_v) < (double) EPS_DOUBLE)
  1627.       continue;
  1628.  
  1629.     temperature_v = current_user_parameter_temp[index_v];
  1630.     parameter_v = last_saved_state->parameter[index_v];
  1631.  
  1632.     /* Handle discrete integer parameters. */
  1633.     if (INTEGER_PARAMETER (index_v))
  1634.       {
  1635.     min_parameter_v -= HALF;
  1636.     max_parameter_v += HALF;
  1637.     parameter_range_v = max_parameter_v - min_parameter_v;
  1638.       }
  1639.     /* generate a new state x within the parameter bounds */
  1640.     for (;;)
  1641.       {
  1642.     x = parameter_v
  1643.       + generate_asa_state (user_random_generator, &temperature_v)
  1644.       * parameter_range_v;
  1645.     /* exit the loop if within its valid parameter range */
  1646.     if (x <= max_parameter_v - (double) EPS_DOUBLE
  1647.         && x >= min_parameter_v + (double) EPS_DOUBLE)
  1648.       break;
  1649.       }
  1650.  
  1651.     /* Handle discrete integer parameters.
  1652.        You might have to check rounding on your machine. */
  1653.     if (INTEGER_PARAMETER (index_v))
  1654.       {
  1655.     if (x < min_parameter_v + HALF)
  1656.       x = min_parameter_v + HALF + (double) EPS_DOUBLE;
  1657.     if (x > max_parameter_v - HALF)
  1658.       x = max_parameter_v - HALF + (double) EPS_DOUBLE;
  1659.  
  1660.     if (x + HALF > ZERO)
  1661.       {
  1662.         x = (int) (x + HALF);
  1663.       }
  1664.     else
  1665.       {
  1666.         x = (int) (x - HALF);
  1667.       }
  1668.     if (x > parameter_maximum[index_v])
  1669.       x = parameter_maximum[index_v];
  1670.     if (x < parameter_minimum[index_v])
  1671.       x = parameter_minimum[index_v];
  1672.       }
  1673.     /* save the newly generated value */
  1674.     current_generated_state->parameter[index_v] = x;
  1675.  
  1676.     if (OPTIONS->SEQUENTIAL_PARAMETERS >= 0)
  1677.       break;
  1678.   }
  1679.  
  1680. }
  1681.  
  1682. /***********************************************************************
  1683. * generate_asa_state
  1684. *       This function generates a single value according to the
  1685. *       ASA generating function and the passed temperature
  1686. ***********************************************************************/
  1687. #if HAVE_ANSI
  1688. double
  1689. generate_asa_state (double (*user_random_generator) (),
  1690.             double *temp)
  1691. #else
  1692. double
  1693. generate_asa_state (user_random_generator,
  1694.             temp)
  1695.      double (*user_random_generator) ();
  1696.      double *temp;
  1697. #endif
  1698. {
  1699.   double x, y, z;
  1700.  
  1701.   x = (*user_random_generator) ();
  1702.   y = x < HALF ? -ONE : ONE;
  1703.   z = y * *temp * (pow ((ONE + ONE / *temp), fabs (TWO * x - ONE)) - ONE);
  1704.  
  1705.   return (z);
  1706.  
  1707. }
  1708.  
  1709. /***********************************************************************
  1710. * accept_new_state
  1711. *    This procedure accepts or rejects a newly generated state,
  1712. *    depending on whether the difference between new and old
  1713. *    cost functions passes a statistical test. If accepted,
  1714. *    the current state is updated.
  1715. ***********************************************************************/
  1716. #if HAVE_ANSI
  1717. void
  1718. accept_new_state (double (*user_random_generator) (),
  1719.           double *parameter_minimum,
  1720.           double *parameter_maximum,
  1721.           double *current_cost_temperature,
  1722. #if ASA_SAMPLE
  1723.           double *current_user_parameter_temp,
  1724. #endif
  1725.           ALLOC_INT * number_parameters,
  1726.           LONG_INT * recent_number_acceptances,
  1727.           LONG_INT * number_accepted,
  1728.           LONG_INT * index_cost_acceptances,
  1729.           LONG_INT * number_acceptances_saved,
  1730.           LONG_INT * recent_number_generated,
  1731.           LONG_INT * number_generated,
  1732.           LONG_INT * index_parameter_generations,
  1733.           STATE * current_generated_state,
  1734.           STATE * last_saved_state,
  1735. #if ASA_SAMPLE
  1736.           FILE * ptr_asa_out,
  1737. #endif
  1738.           USER_DEFINES * OPTIONS)
  1739. #else
  1740. void
  1741. accept_new_state (user_random_generator,
  1742.           parameter_minimum,
  1743.           parameter_maximum,
  1744.           current_cost_temperature,
  1745. #if ASA_SAMPLE
  1746.           current_user_parameter_temp,
  1747. #endif
  1748.           number_parameters,
  1749.           recent_number_acceptances,
  1750.           number_accepted,
  1751.           index_cost_acceptances,
  1752.           number_acceptances_saved,
  1753.           recent_number_generated,
  1754.           number_generated,
  1755.           index_parameter_generations,
  1756.           current_generated_state,
  1757.           last_saved_state,
  1758. #if ASA_SAMPLE
  1759.           ptr_asa_out,
  1760. #endif
  1761.           OPTIONS)
  1762.      double (*user_random_generator) ();
  1763.      double *parameter_minimum;
  1764.      double *parameter_maximum;
  1765.      double *current_cost_temperature;
  1766. #if ASA_SAMPLE
  1767.      double *current_user_parameter_temp;
  1768. #endif
  1769.      ALLOC_INT *number_parameters;
  1770.      LONG_INT *recent_number_acceptances;
  1771.      LONG_INT *number_accepted;
  1772.      LONG_INT *index_cost_acceptances;
  1773.      LONG_INT *number_acceptances_saved;
  1774.      LONG_INT *recent_number_generated;
  1775.      LONG_INT *number_generated;
  1776.      LONG_INT *index_parameter_generations;
  1777.      STATE *current_generated_state;
  1778.      STATE *last_saved_state;
  1779. #if ASA_SAMPLE
  1780.      FILE *ptr_asa_out;
  1781. #endif
  1782.      USER_DEFINES *OPTIONS;
  1783.  
  1784. #endif
  1785. {
  1786.   double delta_cost, boltz_test, unif_test;
  1787.   ALLOC_INT index_v;
  1788.   double curr_cost_temp;
  1789. #if ASA_SAMPLE
  1790.   LONG_INT active_params;
  1791.   double weight_param_ind, weight_aver, range;
  1792. #endif
  1793.  
  1794.   /* update accepted and generated count */
  1795.   ++*number_acceptances_saved;
  1796.   ++*recent_number_generated;
  1797.   ++*number_generated;
  1798.  
  1799.   /* increment the parameter index generation for each parameter */
  1800.   if (OPTIONS->SEQUENTIAL_PARAMETERS >= 0)
  1801.     {
  1802.       /* ignore parameters with too small a range */
  1803.       if (!PARAMETER_RANGE_TOO_SMALL (OPTIONS->SEQUENTIAL_PARAMETERS))
  1804.     ++index_parameter_generations[OPTIONS->SEQUENTIAL_PARAMETERS];
  1805.     }
  1806.   else
  1807.     {
  1808.       VFOR (index_v)
  1809.       {
  1810.     if (!PARAMETER_RANGE_TOO_SMALL (index_v))
  1811.       ++index_parameter_generations[index_v];
  1812.       }
  1813.     }
  1814.  
  1815.   /* effective cost function for testing acceptance criteria,
  1816.      calculate the cost difference and divide by the temperature */
  1817. #if USER_COST_SCHEDULE
  1818.   curr_cost_temp =
  1819.     (OPTIONS->cost_schedule (*current_cost_temperature, OPTIONS)
  1820.      + EPS_DOUBLE);
  1821. #else
  1822.   curr_cost_temp = *current_cost_temperature;
  1823. #endif
  1824.   delta_cost = (current_generated_state->cost - last_saved_state->cost)
  1825.     / (curr_cost_temp + EPS_DOUBLE);
  1826.   boltz_test = MIN (ONE, (exp (EXPONENT_CHECK (-delta_cost))));
  1827.   unif_test = (*user_random_generator) ();
  1828.  
  1829. #if ASA_SAMPLE
  1830.   active_params = 0;
  1831.   weight_aver = ZERO;
  1832.   VFOR (index_v)
  1833.   {
  1834.     /* ignore parameters with too small a range */
  1835.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1836.       continue;
  1837.     ++active_params;
  1838.     range = parameter_maximum[index_v] - parameter_minimum[index_v];
  1839.     weight_param_ind = TWO
  1840.       * (fabs ((last_saved_state->parameter[index_v]
  1841.          - current_generated_state->parameter[index_v]) / (TWO * range))
  1842.      + current_user_parameter_temp[index_v])
  1843.       * log (ONE + ONE / current_user_parameter_temp[index_v]);
  1844.     weight_aver += weight_param_ind;
  1845.     OPTIONS->bias_generated[index_v] = ONE / weight_param_ind;
  1846.   }
  1847.   weight_aver /= (double) active_params;
  1848.   OPTIONS->average_weights = weight_aver;
  1849.   OPTIONS->bias_acceptance = boltz_test;
  1850.   if (boltz_test >= unif_test)
  1851.     {
  1852.       ++OPTIONS->n_accepted;
  1853.     }
  1854.  
  1855. #if ASA_PRINT
  1856.   if (OPTIONS->limit_weights < OPTIONS->average_weights)
  1857.     {
  1858.       fprintf (ptr_asa_out, ":SAMPLE#\n");
  1859.       if (boltz_test >= unif_test)
  1860.     {
  1861. #if INT_LONG
  1862.       fprintf (ptr_asa_out, ":SAMPLE+ %10ld %12.7g %12.7g %12.7g %12.7g\n",
  1863. #else
  1864.       fprintf (ptr_asa_out, ":SAMPLE+ %10d %12.7g %12.7g %12.7g\n",
  1865. #endif
  1866.            OPTIONS->n_accepted, current_generated_state->cost,
  1867.            curr_cost_temp, OPTIONS->bias_acceptance,
  1868.            OPTIONS->average_weights);
  1869.       VFOR (index_v)
  1870.       {
  1871.         /* ignore parameters with too small a range */
  1872.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1873.           continue;
  1874.         range = parameter_maximum[index_v] - parameter_minimum[index_v];
  1875. #if INT_ALLOC
  1876.         fprintf (ptr_asa_out, ":SAMPLE %11d %12.7g %12.7g %12.7g %12.7g\n",
  1877. #else
  1878. #if INT_LONG
  1879.         fprintf (ptr_asa_out, ":SAMPLE %11ld %12.7g %12.7g %12.7g %12.7g\n",
  1880. #else
  1881.         fprintf (ptr_asa_out, ":SAMPLE %11d %12.7g %12.7g %12.7g %12.7g\n",
  1882. #endif
  1883. #endif
  1884.              index_v, current_generated_state->parameter[index_v],
  1885.              current_user_parameter_temp[index_v],
  1886.              OPTIONS->bias_generated[index_v], range);
  1887.       }
  1888.     }
  1889.       else
  1890.     {
  1891. #if INT_LONG
  1892.       fprintf (ptr_asa_out, ":SAMPLE %11ld %12.7g %12.7g %12.7g %12.7g\n",
  1893. #else
  1894.       fprintf (ptr_asa_out, ":SAMPLE %11d %12.7g %12.7g %12.7g\n",
  1895. #endif
  1896.            OPTIONS->n_accepted, last_saved_state->cost,
  1897.            curr_cost_temp, OPTIONS->bias_acceptance,
  1898.            OPTIONS->average_weights);
  1899.       VFOR (index_v)
  1900.       {
  1901.         /* ignore parameters with too small a range */
  1902.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1903.           continue;
  1904.         range = parameter_maximum[index_v] - parameter_minimum[index_v];
  1905. #if INT_ALLOC
  1906.         fprintf (ptr_asa_out, ":SAMPLE %11d %12.7g %12.7g %12.7g %12.7g\n",
  1907. #else
  1908. #if INT_LONG
  1909.         fprintf (ptr_asa_out, ":SAMPLE %11ld %12.7g %12.7g %12.7g %12.7g\n",
  1910. #else
  1911.         fprintf (ptr_asa_out, ":SAMPLE %11d %12.7g %12.7g %12.7g %12.7g\n",
  1912. #endif
  1913. #endif
  1914.              index_v, last_saved_state->parameter[index_v],
  1915.              current_user_parameter_temp[index_v],
  1916.              OPTIONS->bias_generated[index_v], range);
  1917.       }
  1918.     }
  1919.     }
  1920. #endif
  1921. #endif /* ASA_SAMPLE */
  1922.  
  1923.   /* accept/reject the new state */
  1924.   if (boltz_test >= unif_test)
  1925.     {
  1926.       /* copy current state to the last saved state */
  1927.  
  1928.       last_saved_state->cost = current_generated_state->cost;
  1929.       VFOR (index_v)
  1930.       {
  1931.     /* ignore parameters with too small a range */
  1932.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  1933.       continue;
  1934.     last_saved_state->parameter[index_v] =
  1935.       current_generated_state->parameter[index_v];
  1936.       }
  1937.  
  1938.       /* update acceptance counts */
  1939.       ++*recent_number_acceptances;
  1940.       ++*number_accepted;
  1941.       ++*index_cost_acceptances;
  1942.       *number_acceptances_saved = *number_accepted;
  1943.     }
  1944. }
  1945.  
  1946. /***********************************************************************
  1947. * reanneal
  1948. *    Readjust temperatures of generating and acceptance functions
  1949. ***********************************************************************/
  1950. #if HAVE_ANSI
  1951. void
  1952. reanneal (double *parameter_minimum,
  1953.       double *parameter_maximum,
  1954.       double *tangents,
  1955.       double *maximum_tangent,
  1956.       double *current_cost_temperature,
  1957.       double *initial_cost_temperature,
  1958.       double *temperature_scale_cost,
  1959.       double *current_user_parameter_temp,
  1960.       double *initial_user_parameter_temp,
  1961.       double *temperature_scale_parameters,
  1962.       double *quench_param,
  1963.       double *quench_cost,
  1964.       ALLOC_INT * number_parameters,
  1965.       int *parameter_type,
  1966.       LONG_INT * index_cost_acceptances,
  1967.       LONG_INT * index_parameter_generations,
  1968.       STATE * last_saved_state,
  1969.       STATE * best_generated_state,
  1970.       USER_DEFINES * OPTIONS)
  1971. #else
  1972. void
  1973. reanneal (parameter_minimum,
  1974.       parameter_maximum,
  1975.       tangents,
  1976.       maximum_tangent,
  1977.       current_cost_temperature,
  1978.       initial_cost_temperature,
  1979.       temperature_scale_cost,
  1980.       current_user_parameter_temp,
  1981.       initial_user_parameter_temp,
  1982.       temperature_scale_parameters,
  1983.       quench_param,
  1984.       quench_cost,
  1985.       number_parameters,
  1986.       parameter_type,
  1987.       index_cost_acceptances,
  1988.       index_parameter_generations,
  1989.       last_saved_state,
  1990.       best_generated_state,
  1991.       OPTIONS)
  1992.      double *parameter_minimum;
  1993.      double *parameter_maximum;
  1994.      double *tangents;
  1995.      double *maximum_tangent;
  1996.      double *current_cost_temperature;
  1997.      double *initial_cost_temperature;
  1998.      double *temperature_scale_cost;
  1999.      double *current_user_parameter_temp;
  2000.      double *initial_user_parameter_temp;
  2001.      double *temperature_scale_parameters;
  2002.      double *quench_param;
  2003.      double *quench_cost;
  2004.      ALLOC_INT *number_parameters;
  2005.      int *parameter_type;
  2006.      LONG_INT *index_cost_acceptances;
  2007.      LONG_INT *index_parameter_generations;
  2008.      STATE *last_saved_state;
  2009.      STATE *best_generated_state;
  2010.      USER_DEFINES *OPTIONS;
  2011. #endif
  2012. {
  2013.   ALLOC_INT index_v;
  2014.   double tmp_var_db3;
  2015.   double new_temperature;    /* the temperature */
  2016.   double log_new_temperature_ratio;
  2017.   double log_init_cur_temp_ratio;
  2018.   double temperature_rescale_power;
  2019.   double tmp_dbl;
  2020.  
  2021.   VFOR (index_v)
  2022.   {
  2023.     if (NO_REANNEAL (index_v))
  2024.       continue;
  2025.  
  2026.     /* use the temp double to prevent overflow */
  2027.     tmp_dbl = (double) index_parameter_generations[index_v];
  2028.  
  2029.     /* skip parameters with too small range or integer parameters */
  2030.     if (OPTIONS->INCLUDE_INTEGER_PARAMETERS == TRUE)
  2031.       {
  2032.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2033.       continue;
  2034.       }
  2035.     else
  2036.       {
  2037.     if (PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2038.         INTEGER_PARAMETER (index_v))
  2039.       continue;
  2040.       }
  2041.  
  2042.     /* ignore parameters with too small tangents */
  2043.     if (fabs (tangents[index_v]) < (double) EPS_DOUBLE)
  2044.       continue;
  2045.  
  2046.     /* reset the index of parameter generations appropriately */
  2047. #if USER_REANNEAL_FUNCTION
  2048.     new_temperature =
  2049.       fabs (OPTIONS->reanneal_function (current_user_parameter_temp[index_v],
  2050.                     tangents[index_v],
  2051.                     *maximum_tangent,
  2052.                     OPTIONS));
  2053. #else
  2054.     new_temperature =
  2055.       fabs (REANNEAL_FUNCTION (current_user_parameter_temp[index_v],
  2056.                    tangents[index_v],
  2057.                    *maximum_tangent));
  2058. #endif
  2059.     if (new_temperature < initial_user_parameter_temp[index_v])
  2060.       {
  2061.     log_init_cur_temp_ratio =
  2062.       fabs (log (((double) EPS_DOUBLE
  2063.               + initial_user_parameter_temp[index_v])
  2064.              / ((double) EPS_DOUBLE + new_temperature)));
  2065.     tmp_dbl = (double) EPS_DOUBLE
  2066.       + pow (log_init_cur_temp_ratio
  2067.          / temperature_scale_parameters[index_v],
  2068.          (double) *number_parameters
  2069.          / quench_param[index_v]);
  2070.       }
  2071.     else
  2072.       {
  2073.     tmp_dbl = ONE;
  2074.       }
  2075.  
  2076.     /* Reset index_parameter_generations if index reset too large,
  2077.        and also reset the initial_user_parameter_temp, to achieve
  2078.        the same new temperature. */
  2079.     while (tmp_dbl
  2080.        > ((double) OPTIONS->MAXIMUM_REANNEAL_INDEX))
  2081.       {
  2082.     log_new_temperature_ratio =
  2083.       -temperature_scale_parameters[index_v] *
  2084.       pow (tmp_dbl,
  2085.            quench_param[index_v] / (double) *number_parameters);
  2086.     log_new_temperature_ratio =
  2087.       EXPONENT_CHECK (log_new_temperature_ratio);
  2088.     new_temperature =
  2089.       initial_user_parameter_temp[index_v]
  2090.       * exp (log_new_temperature_ratio);
  2091.     tmp_dbl /=
  2092.       OPTIONS->REANNEAL_RESCALE;
  2093.     temperature_rescale_power = ONE
  2094.       / pow (OPTIONS->REANNEAL_RESCALE,
  2095.          quench_param[index_v]
  2096.          / (double) *number_parameters);
  2097.     initial_user_parameter_temp[index_v] =
  2098.       new_temperature * pow (initial_user_parameter_temp
  2099.                  [index_v] / new_temperature,
  2100.                  temperature_rescale_power);
  2101.       }
  2102.     /* restore from temporary double */
  2103.     index_parameter_generations[index_v] = (LONG_INT) tmp_dbl;
  2104.   }
  2105.  
  2106.   /* reanneal : Reset the index of cost acceptances to take
  2107.      into account finer detail in cost terrain. */
  2108.  
  2109.   /* set the starting cost_temperature to last cost found so far */
  2110.   if (*initial_cost_temperature > fabs (last_saved_state->cost))
  2111.     *initial_cost_temperature = fabs (last_saved_state->cost);
  2112.  
  2113.   tmp_dbl = (double) *index_cost_acceptances;
  2114.  
  2115.   if (*current_cost_temperature > fabs (best_generated_state->cost))
  2116.     {
  2117.       tmp_var_db3 =
  2118.     fabs (log (((double) EPS_DOUBLE + *initial_cost_temperature) /
  2119.         ((double) EPS_DOUBLE + fabs (best_generated_state->cost))));
  2120.       tmp_dbl = (double) EPS_DOUBLE
  2121.     + pow (tmp_var_db3
  2122.            / *temperature_scale_cost,
  2123.            (double) *number_parameters / (*quench_cost));
  2124.     }
  2125.   else
  2126.     {
  2127.       log_init_cur_temp_ratio =
  2128.     fabs (log (((double) EPS_DOUBLE + *initial_cost_temperature) /
  2129.            ((double) EPS_DOUBLE + *current_cost_temperature)));
  2130.       tmp_dbl = (double) EPS_DOUBLE
  2131.     + pow (log_init_cur_temp_ratio
  2132.            / *temperature_scale_cost,
  2133.            (double) *number_parameters / (*quench_cost));
  2134.     }
  2135.  
  2136.   /* reset index_cost_temperature if index reset too large */
  2137.   while (tmp_dbl > ((double) OPTIONS->MAXIMUM_REANNEAL_INDEX))
  2138.     {
  2139.       log_new_temperature_ratio = -*temperature_scale_cost *
  2140.     pow ((double) tmp_dbl,
  2141.          (*quench_cost) / (double) *number_parameters);
  2142.       log_new_temperature_ratio =
  2143.     EXPONENT_CHECK (log_new_temperature_ratio);
  2144.       new_temperature = *initial_cost_temperature
  2145.     * exp (log_new_temperature_ratio);
  2146.       tmp_dbl /= OPTIONS->REANNEAL_RESCALE;
  2147.       temperature_rescale_power = ONE
  2148.     / pow (OPTIONS->REANNEAL_RESCALE,
  2149.            (*quench_cost)
  2150.            / (double) *number_parameters);
  2151.       *initial_cost_temperature =
  2152.     new_temperature * pow (*initial_cost_temperature
  2153.                    / new_temperature,
  2154.                    temperature_rescale_power);
  2155.     }
  2156.   *index_cost_acceptances = (LONG_INT) tmp_dbl;
  2157. }
  2158.  
  2159. /***********************************************************************
  2160. * cost_derivatives
  2161. *    This procedure calculates the derivatives of the cost function
  2162. *    with respect to its parameters.  The first derivatives are
  2163. *    used as a sensitivity measure for reannealing.  The second
  2164. *    derivatives are calculated only if *curvature_flag=TRUE;
  2165. *    these are a measure of the covariance of the fit when a
  2166. *    minimum is found.
  2167. ***********************************************************************/
  2168.       /* Calculate the numerical derivatives of the best
  2169.          generated state found so far */
  2170.  
  2171.       /* In this implementation of ASA, no checks are made for
  2172.          *valid_state_generated_flag=FALSE for differential neighbors
  2173.          to the current best state. */
  2174.  
  2175.       /* Assuming no information is given about the metric of the parameter
  2176.          space, use simple Cartesian space to calculate curvatures. */
  2177.  
  2178. #if HAVE_ANSI
  2179. void
  2180. cost_derivatives (double (*user_cost_function) (
  2181.                double *, double *, double *, double *, double *,
  2182.               ALLOC_INT *, int *, int *, int *, USER_DEFINES *),
  2183.           double *parameter_minimum,
  2184.           double *parameter_maximum,
  2185.           double *tangents,
  2186.           double *curvature,
  2187.           double *maximum_tangent,
  2188.           double *delta_parameter,
  2189.           ALLOC_INT * number_parameters,
  2190.           int *parameter_type,
  2191.           int *exit_status,
  2192.           int *curvature_flag,
  2193.           int *valid_state_generated_flag,
  2194.           LONG_INT * number_invalid_generated_states,
  2195.           STATE * current_generated_state,
  2196.           STATE * best_generated_state,
  2197.           FILE * ptr_asa_out,
  2198.           USER_DEFINES * OPTIONS)
  2199. #else
  2200. void
  2201. cost_derivatives (user_cost_function,
  2202.           parameter_minimum,
  2203.           parameter_maximum,
  2204.           tangents,
  2205.           curvature,
  2206.           maximum_tangent,
  2207.           delta_parameter,
  2208.           number_parameters,
  2209.           parameter_type,
  2210.           exit_status,
  2211.           curvature_flag,
  2212.           valid_state_generated_flag,
  2213.           number_invalid_generated_states,
  2214.           current_generated_state,
  2215.           best_generated_state,
  2216.           ptr_asa_out,
  2217.           OPTIONS)
  2218.      double (*user_cost_function) ();
  2219.      double *parameter_minimum;
  2220.      double *parameter_maximum;
  2221.      double *tangents;
  2222.      double *curvature;
  2223.      double *maximum_tangent;
  2224.      double *delta_parameter;
  2225.      ALLOC_INT *number_parameters;
  2226.      int *parameter_type;
  2227.      int *exit_status;
  2228.      int *curvature_flag;
  2229.      int *valid_state_generated_flag;
  2230.      LONG_INT *number_invalid_generated_states;
  2231.      STATE *current_generated_state;
  2232.      STATE *best_generated_state;
  2233.      FILE *ptr_asa_out;
  2234.      USER_DEFINES *OPTIONS;
  2235. #endif
  2236. {
  2237.   ALLOC_INT index_v, index_vv, index_v_vv, index_vv_v;
  2238.   LONG_INT saved_num_invalid_gen_states, tmp_saved;
  2239.   double parameter_v, parameter_vv, parameter_v_offset, parameter_vv_offset;
  2240.   double recent_best_cost;
  2241.   double new_cost_state_1, new_cost_state_2, new_cost_state_3;
  2242.   double delta_parameter_v, delta_parameter_vv;
  2243.  
  2244.   if (OPTIONS->CURVATURE_0 == TRUE)
  2245.     *curvature_flag = FALSE;
  2246.   if (OPTIONS->CURVATURE_0 == -1)
  2247.     *curvature_flag = TRUE;
  2248.  
  2249.   /* save the best cost */
  2250.   recent_best_cost = best_generated_state->cost;
  2251.  
  2252.   /* copy the best state into the current state */
  2253.   VFOR (index_v)
  2254.   {
  2255.     /* ignore parameters with too small ranges */
  2256.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2257.       continue;
  2258.     current_generated_state->parameter[index_v] =
  2259.       best_generated_state->parameter[index_v];
  2260.   }
  2261.  
  2262.   saved_num_invalid_gen_states = (*number_invalid_generated_states);
  2263.  
  2264.   /* set parameters (& possibly constraints) to best state */
  2265.   *valid_state_generated_flag = TRUE;
  2266.   current_generated_state->cost =
  2267.     user_cost_function (current_generated_state->parameter,
  2268.             parameter_minimum,
  2269.             parameter_maximum,
  2270.             tangents,
  2271.             curvature,
  2272.             number_parameters,
  2273.             parameter_type,
  2274.             valid_state_generated_flag,
  2275.             exit_status,
  2276.             OPTIONS);
  2277.   if (*valid_state_generated_flag == FALSE)
  2278.     ++(*number_invalid_generated_states);
  2279.  
  2280.   if (OPTIONS->USER_TANGENTS == TRUE)
  2281.     {
  2282.       *valid_state_generated_flag = FALSE;
  2283.       current_generated_state->cost =
  2284.     user_cost_function (current_generated_state->parameter,
  2285.                 parameter_minimum,
  2286.                 parameter_maximum,
  2287.                 tangents,
  2288.                 curvature,
  2289.                 number_parameters,
  2290.                 parameter_type,
  2291.                 valid_state_generated_flag,
  2292.                 exit_status,
  2293.                 OPTIONS);
  2294.       if (*valid_state_generated_flag == FALSE)
  2295.     ++(*number_invalid_generated_states);
  2296.     }
  2297.   else
  2298.     {
  2299.       /* calculate tangents */
  2300.       VFOR (index_v)
  2301.       {
  2302.     if (NO_REANNEAL (index_v))
  2303.       {
  2304.         tangents[index_v] = ZERO;
  2305.         continue;
  2306.       }
  2307.     /* skip parameters with too small range or integer parameters */
  2308.     if (OPTIONS->INCLUDE_INTEGER_PARAMETERS == TRUE)
  2309.       {
  2310.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2311.           {
  2312.         tangents[index_v] = ZERO;
  2313.         continue;
  2314.           }
  2315.       }
  2316.     else
  2317.       {
  2318.         if (PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2319.         INTEGER_PARAMETER (index_v))
  2320.           {
  2321.         tangents[index_v] = ZERO;
  2322.         continue;
  2323.           }
  2324.       }
  2325.  
  2326.     /* save the v_th parameter and delta_parameter */
  2327.     parameter_v = best_generated_state->parameter[index_v];
  2328.     delta_parameter_v = delta_parameter[index_v];
  2329.  
  2330.     parameter_v_offset = (ONE + delta_parameter_v) * parameter_v;
  2331.     if (parameter_v_offset > parameter_maximum[index_v] ||
  2332.         parameter_v_offset < parameter_minimum[index_v])
  2333.       {
  2334.         delta_parameter_v = -delta_parameter_v;
  2335.         parameter_v_offset = (ONE + delta_parameter_v) * parameter_v;
  2336.       }
  2337.  
  2338.     /* generate the first sample point */
  2339.     current_generated_state->parameter[index_v] = parameter_v_offset;
  2340.     *valid_state_generated_flag = TRUE;
  2341.     current_generated_state->cost =
  2342.       user_cost_function (current_generated_state->parameter,
  2343.                   parameter_minimum,
  2344.                   parameter_maximum,
  2345.                   tangents,
  2346.                   curvature,
  2347.                   number_parameters,
  2348.                   parameter_type,
  2349.                   valid_state_generated_flag,
  2350.                   exit_status,
  2351.                   OPTIONS);
  2352.     if (*valid_state_generated_flag == FALSE)
  2353.       ++(*number_invalid_generated_states);
  2354.     new_cost_state_1 = current_generated_state->cost;
  2355.  
  2356.     /* restore the parameter state */
  2357.     current_generated_state->parameter[index_v] = parameter_v;
  2358.  
  2359.     /* calculate the numerical derivative */
  2360.     tangents[index_v] = (new_cost_state_1 - recent_best_cost)
  2361.       / (delta_parameter_v * parameter_v + (double) EPS_DOUBLE);
  2362.  
  2363.       }
  2364.     }
  2365.  
  2366.   /* find the maximum |tangent| from all tangents */
  2367.   *maximum_tangent = 0;
  2368.   VFOR (index_v)
  2369.   {
  2370.     if (NO_REANNEAL (index_v))
  2371.       continue;
  2372.  
  2373.     /* ignore too small ranges and integer parameters types */
  2374.     if (OPTIONS->INCLUDE_INTEGER_PARAMETERS == TRUE)
  2375.       {
  2376.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2377.       continue;
  2378.       }
  2379.     else
  2380.       {
  2381.     if (PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2382.         INTEGER_PARAMETER (index_v))
  2383.       continue;
  2384.       }
  2385.  
  2386.     /* find the maximum |tangent| (from all tangents) */
  2387.     if (fabs (tangents[index_v]) > *maximum_tangent)
  2388.       {
  2389.     *maximum_tangent = fabs (tangents[index_v]);
  2390.       }
  2391.   }
  2392.  
  2393.   if (*curvature_flag == TRUE || *curvature_flag == -1)
  2394.     {
  2395.       /* calculate diagonal curvatures */
  2396.       VFOR (index_v)
  2397.       {
  2398.     if (NO_REANNEAL (index_v))
  2399.       {
  2400.         index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2401.         curvature[index_v_vv] = ZERO;
  2402.         continue;
  2403.       }
  2404.     /* skip parameters with too small range or integer parameters */
  2405.     if (OPTIONS->INCLUDE_INTEGER_PARAMETERS == TRUE)
  2406.       {
  2407.         if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2408.           {
  2409.         index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2410.         curvature[index_v_vv] = ZERO;
  2411.         continue;
  2412.           }
  2413.         else
  2414.           {
  2415.         if (PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2416.             INTEGER_PARAMETER (index_v))
  2417.           {
  2418.             index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2419.             curvature[index_v_vv] = ZERO;
  2420.             continue;
  2421.           }
  2422.           }
  2423.       }
  2424.  
  2425.     /* save the v_th parameter and delta_parameter */
  2426.     parameter_v = best_generated_state->parameter[index_v];
  2427.     delta_parameter_v = delta_parameter[index_v];
  2428.  
  2429.     if (parameter_v + delta_parameter_v * fabs (parameter_v)
  2430.         > parameter_maximum[index_v])
  2431.       {
  2432.         /* generate the first sample point */
  2433.         current_generated_state->parameter[index_v] =
  2434.           parameter_v - TWO * delta_parameter_v * fabs (parameter_v);
  2435.         *valid_state_generated_flag = TRUE;
  2436.         current_generated_state->cost =
  2437.           user_cost_function (current_generated_state->parameter,
  2438.                   parameter_minimum,
  2439.                   parameter_maximum,
  2440.                   tangents,
  2441.                   curvature,
  2442.                   number_parameters,
  2443.                   parameter_type,
  2444.                   valid_state_generated_flag,
  2445.                   exit_status,
  2446.                   OPTIONS);
  2447.         if (*valid_state_generated_flag == FALSE)
  2448.           ++(*number_invalid_generated_states);
  2449.         new_cost_state_1 = current_generated_state->cost;
  2450.  
  2451.         /* generate the second sample point */
  2452.         current_generated_state->parameter[index_v] =
  2453.           parameter_v - delta_parameter_v * fabs (parameter_v);
  2454.         *valid_state_generated_flag = TRUE;
  2455.         current_generated_state->cost =
  2456.           user_cost_function (current_generated_state->parameter,
  2457.                   parameter_minimum,
  2458.                   parameter_maximum,
  2459.                   tangents,
  2460.                   curvature,
  2461.                   number_parameters,
  2462.                   parameter_type,
  2463.                   valid_state_generated_flag,
  2464.                   exit_status,
  2465.                   OPTIONS);
  2466.         if (*valid_state_generated_flag == FALSE)
  2467.           ++(*number_invalid_generated_states);
  2468.         new_cost_state_2 = current_generated_state->cost;
  2469.  
  2470.         /* restore the parameter state */
  2471.         current_generated_state->parameter[index_v] =
  2472.           parameter_v;
  2473.  
  2474.         /* index_v_vv: row index_v, column index_v */
  2475.         index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2476.  
  2477.         /* calculate and store the curvature */
  2478.         curvature[index_v_vv] =
  2479.           (recent_best_cost - TWO * new_cost_state_2
  2480.            + new_cost_state_1) / (delta_parameter_v * delta_parameter_v
  2481.              * parameter_v * parameter_v + (double) EPS_DOUBLE);
  2482.       }
  2483.     else if (parameter_v - delta_parameter_v * fabs (parameter_v)
  2484.          < parameter_minimum[index_v])
  2485.       {
  2486.         /* generate the first sample point */
  2487.         current_generated_state->parameter[index_v] =
  2488.           parameter_v + TWO * delta_parameter_v * fabs (parameter_v);
  2489.         *valid_state_generated_flag = TRUE;
  2490.         current_generated_state->cost =
  2491.           user_cost_function (current_generated_state->parameter,
  2492.                   parameter_minimum,
  2493.                   parameter_maximum,
  2494.                   tangents,
  2495.                   curvature,
  2496.                   number_parameters,
  2497.                   parameter_type,
  2498.                   valid_state_generated_flag,
  2499.                   exit_status,
  2500.                   OPTIONS);
  2501.         if (*valid_state_generated_flag == FALSE)
  2502.           ++(*number_invalid_generated_states);
  2503.         new_cost_state_1 = current_generated_state->cost;
  2504.  
  2505.         /* generate the second sample point */
  2506.         current_generated_state->parameter[index_v] =
  2507.           parameter_v + delta_parameter_v * fabs (parameter_v);
  2508.         *valid_state_generated_flag = TRUE;
  2509.         current_generated_state->cost =
  2510.           user_cost_function (current_generated_state->parameter,
  2511.                   parameter_minimum,
  2512.                   parameter_maximum,
  2513.                   tangents,
  2514.                   curvature,
  2515.                   number_parameters,
  2516.                   parameter_type,
  2517.                   valid_state_generated_flag,
  2518.                   exit_status,
  2519.                   OPTIONS);
  2520.         if (*valid_state_generated_flag == FALSE)
  2521.           ++(*number_invalid_generated_states);
  2522.         new_cost_state_2 = current_generated_state->cost;
  2523.  
  2524.         /* restore the parameter state */
  2525.         current_generated_state->parameter[index_v] =
  2526.           parameter_v;
  2527.  
  2528.         /* index_v_vv: row index_v, column index_v */
  2529.         index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2530.  
  2531.         /* calculate and store the curvature */
  2532.         curvature[index_v_vv] =
  2533.           (recent_best_cost - TWO * new_cost_state_2
  2534.            + new_cost_state_1) / (delta_parameter_v * delta_parameter_v
  2535.              * parameter_v * parameter_v + (double) EPS_DOUBLE);
  2536.       }
  2537.     else
  2538.       {
  2539.         /* generate the first sample point */
  2540.         parameter_v_offset = (ONE + delta_parameter_v) * parameter_v;
  2541.         current_generated_state->parameter[index_v] = parameter_v_offset;
  2542.         *valid_state_generated_flag = TRUE;
  2543.         current_generated_state->cost =
  2544.           user_cost_function (current_generated_state->parameter,
  2545.                   parameter_minimum,
  2546.                   parameter_maximum,
  2547.                   tangents,
  2548.                   curvature,
  2549.                   number_parameters,
  2550.                   parameter_type,
  2551.                   valid_state_generated_flag,
  2552.                   exit_status,
  2553.                   OPTIONS);
  2554.         if (*valid_state_generated_flag == FALSE)
  2555.           ++(*number_invalid_generated_states);
  2556.         new_cost_state_1 = current_generated_state->cost;
  2557.  
  2558.         /* generate the second sample point */
  2559.         current_generated_state->parameter[index_v] =
  2560.           (ONE - delta_parameter_v) * parameter_v;
  2561.         *valid_state_generated_flag = TRUE;
  2562.         current_generated_state->cost =
  2563.           user_cost_function (current_generated_state->parameter,
  2564.                   parameter_minimum,
  2565.                   parameter_maximum,
  2566.                   tangents,
  2567.                   curvature,
  2568.                   number_parameters,
  2569.                   parameter_type,
  2570.                   valid_state_generated_flag,
  2571.                   exit_status,
  2572.                   OPTIONS);
  2573.         if (*valid_state_generated_flag == FALSE)
  2574.           ++(*number_invalid_generated_states);
  2575.         new_cost_state_2 = current_generated_state->cost;
  2576.  
  2577.         /* restore the parameter state */
  2578.         current_generated_state->parameter[index_v] =
  2579.           parameter_v;
  2580.  
  2581.         /* index_v_vv: row index_v, column index_v */
  2582.         index_v_vv = ROW_COL_INDEX (index_v, index_v);
  2583.  
  2584.         /* calculate and store the curvature */
  2585.         curvature[index_v_vv] =
  2586.           (new_cost_state_2 - TWO * recent_best_cost
  2587.            + new_cost_state_1) / (delta_parameter_v * delta_parameter_v
  2588.              * parameter_v * parameter_v + (double) EPS_DOUBLE);
  2589.       }
  2590.       }
  2591.  
  2592.       /* calculate off-diagonal curvatures */
  2593.       VFOR (index_v)
  2594.       {
  2595.     /* save the v_th parameter and delta_x */
  2596.     parameter_v = current_generated_state->parameter[index_v];
  2597.     delta_parameter_v = delta_parameter[index_v];
  2598.  
  2599.     VFOR (index_vv)
  2600.     {
  2601.       /* index_v_vv: row index_v, column index_vv */
  2602.       index_v_vv = ROW_COL_INDEX (index_v, index_vv);
  2603.       index_vv_v = ROW_COL_INDEX (index_vv, index_v);
  2604.  
  2605.       if (NO_REANNEAL (index_vv) || NO_REANNEAL (index_v))
  2606.         {
  2607.           curvature[index_vv_v] = curvature[index_v_vv] = ZERO;
  2608.           continue;
  2609.         }
  2610.  
  2611.       /* calculate only the upper diagonal */
  2612.       if (index_v <= index_vv)
  2613.         continue;
  2614.  
  2615.       /* skip parms with too small range or integer parameters */
  2616.       if (OPTIONS->INCLUDE_INTEGER_PARAMETERS == TRUE)
  2617.         {
  2618.           if (PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2619.           PARAMETER_RANGE_TOO_SMALL (index_vv))
  2620.         {
  2621.           curvature[index_vv_v] = curvature[index_v_vv] = ZERO;
  2622.           continue;
  2623.         }
  2624.         }
  2625.       else
  2626.         {
  2627.           if (INTEGER_PARAMETER (index_v) ||
  2628.           INTEGER_PARAMETER (index_vv) ||
  2629.           PARAMETER_RANGE_TOO_SMALL (index_v) ||
  2630.           PARAMETER_RANGE_TOO_SMALL (index_vv))
  2631.         {
  2632.           curvature[index_vv_v] = curvature[index_v_vv] = ZERO;
  2633.           continue;
  2634.         }
  2635.         }
  2636.       /* save the vv_th parameter and delta_parameter */
  2637.       parameter_vv =
  2638.         current_generated_state->parameter[index_vv];
  2639.       delta_parameter_vv = delta_parameter[index_vv];
  2640.  
  2641.       /* generate first sample point */
  2642.       parameter_v_offset = current_generated_state->parameter[index_v] =
  2643.         (ONE + delta_parameter_v) * parameter_v;
  2644.       parameter_vv_offset = current_generated_state->parameter[index_vv] =
  2645.         (ONE + delta_parameter_vv) * parameter_vv;
  2646.       if (parameter_v_offset > parameter_maximum[index_v] ||
  2647.           parameter_v_offset < parameter_minimum[index_v])
  2648.         {
  2649.           delta_parameter_v = -delta_parameter_v;
  2650.           current_generated_state->parameter[index_v] =
  2651.         (ONE + delta_parameter_v) * parameter_v;
  2652.         }
  2653.       if (parameter_vv_offset > parameter_maximum[index_vv] ||
  2654.           parameter_vv_offset < parameter_minimum[index_vv])
  2655.         {
  2656.           delta_parameter_vv = -delta_parameter_vv;
  2657.           current_generated_state->parameter[index_vv] =
  2658.         (ONE + delta_parameter_vv) * parameter_vv;
  2659.         }
  2660.  
  2661.       *valid_state_generated_flag = TRUE;
  2662.       current_generated_state->cost =
  2663.         user_cost_function (current_generated_state->parameter,
  2664.                 parameter_minimum,
  2665.                 parameter_maximum,
  2666.                 tangents,
  2667.                 curvature,
  2668.                 number_parameters,
  2669.                 parameter_type,
  2670.                 valid_state_generated_flag,
  2671.                 exit_status,
  2672.                 OPTIONS);
  2673.       if (*valid_state_generated_flag == FALSE)
  2674.         ++(*number_invalid_generated_states);
  2675.       new_cost_state_1 = current_generated_state->cost;
  2676.  
  2677.       /* restore the v_th parameter */
  2678.       current_generated_state->parameter[index_v] =
  2679.         parameter_v;
  2680.  
  2681.       /* generate second sample point */
  2682.       *valid_state_generated_flag = TRUE;
  2683.       current_generated_state->cost =
  2684.         user_cost_function (current_generated_state->parameter,
  2685.                 parameter_minimum,
  2686.                 parameter_maximum,
  2687.                 tangents,
  2688.                 curvature,
  2689.                 number_parameters,
  2690.                 parameter_type,
  2691.                 valid_state_generated_flag,
  2692.                 exit_status,
  2693.                 OPTIONS);
  2694.       if (*valid_state_generated_flag == FALSE)
  2695.         ++(*number_invalid_generated_states);
  2696.       new_cost_state_2 = current_generated_state->cost;
  2697.  
  2698.       /* restore the vv_th parameter */
  2699.       current_generated_state->parameter[index_vv] =
  2700.         parameter_vv;
  2701.  
  2702.       /* generate third sample point */
  2703.       current_generated_state->parameter[index_v] =
  2704.         (ONE + delta_parameter_v) * parameter_v;
  2705.       *valid_state_generated_flag = TRUE;
  2706.       current_generated_state->cost =
  2707.         user_cost_function (current_generated_state->parameter,
  2708.                 parameter_minimum,
  2709.                 parameter_maximum,
  2710.                 tangents,
  2711.                 curvature,
  2712.                 number_parameters,
  2713.                 parameter_type,
  2714.                 valid_state_generated_flag,
  2715.                 exit_status,
  2716.                 OPTIONS);
  2717.       if (*valid_state_generated_flag == FALSE)
  2718.         ++(*number_invalid_generated_states);
  2719.       new_cost_state_3 = current_generated_state->cost;
  2720.  
  2721.       /* restore the v_th parameter */
  2722.       current_generated_state->parameter[index_v] =
  2723.         parameter_v;
  2724.  
  2725.       /* calculate and store the curvature */
  2726.       curvature[index_vv_v] = curvature[index_v_vv] =
  2727.         (new_cost_state_1 - new_cost_state_2
  2728.          - new_cost_state_3 + recent_best_cost)
  2729.         / (delta_parameter_v * delta_parameter_vv
  2730.            * parameter_v * parameter_vv
  2731.            + (double) EPS_DOUBLE);
  2732.     }
  2733.       }
  2734.     }
  2735.   /* restore the best cost function value */
  2736.   current_generated_state->cost = recent_best_cost;
  2737. #if ASA_PRINT
  2738.   tmp_saved = *number_invalid_generated_states - saved_num_invalid_gen_states;
  2739.   if (tmp_saved > 0)
  2740. #if INT_LONG
  2741.     fprintf (ptr_asa_out,
  2742.       "Generated %ld invalid states when calculating the derivatives\n",
  2743.          tmp_saved);
  2744. #else
  2745.     fprintf (ptr_asa_out,
  2746.        "Generated %d invalid states when calculating the derivatives\n",
  2747.          tmp_saved);
  2748. #endif
  2749. #endif /* ASA_PRINT */
  2750.   *number_invalid_generated_states = saved_num_invalid_gen_states;
  2751.  
  2752. }
  2753.  
  2754. #if ASA_PRINT
  2755. /***********************************************************************
  2756. * print_state
  2757. *    Prints a description of the current state of the system
  2758. ***********************************************************************/
  2759. #if HAVE_ANSI
  2760. void
  2761. print_state (double *parameter_minimum,
  2762.          double *parameter_maximum,
  2763.          double *tangents,
  2764.          double *curvature,
  2765.          double *current_cost_temperature,
  2766.          double *current_user_parameter_temp,
  2767.          double *accepted_to_generated_ratio,
  2768.          ALLOC_INT * number_parameters,
  2769.          int *curvature_flag,
  2770.          LONG_INT * number_accepted,
  2771.          LONG_INT * index_cost_acceptances,
  2772.          LONG_INT * number_generated,
  2773.          LONG_INT * number_invalid_generated_states,
  2774.          STATE * last_saved_state,
  2775.          STATE * best_generated_state,
  2776.          FILE * ptr_asa_out,
  2777.          USER_DEFINES * OPTIONS)
  2778. #else
  2779. void
  2780. print_state (parameter_minimum,
  2781.          parameter_maximum,
  2782.          tangents,
  2783.          curvature,
  2784.          current_cost_temperature,
  2785.          current_user_parameter_temp,
  2786.          accepted_to_generated_ratio,
  2787.          number_parameters,
  2788.          curvature_flag,
  2789.          number_accepted,
  2790.          index_cost_acceptances,
  2791.          number_generated,
  2792.          number_invalid_generated_states,
  2793.          last_saved_state,
  2794.          best_generated_state,
  2795.          ptr_asa_out,
  2796.          OPTIONS)
  2797.      double *parameter_minimum;
  2798.      double *parameter_maximum;
  2799.      double *tangents;
  2800.      double *curvature;
  2801.      double *current_cost_temperature;
  2802.      double *current_user_parameter_temp;
  2803.      double *accepted_to_generated_ratio;
  2804.      ALLOC_INT *number_parameters;
  2805.      int *curvature_flag;
  2806.      LONG_INT *number_accepted;
  2807.      LONG_INT *index_cost_acceptances;
  2808.      LONG_INT *number_generated;
  2809.      LONG_INT *number_invalid_generated_states;
  2810.      STATE *last_saved_state;
  2811.      STATE *best_generated_state;
  2812.      FILE *ptr_asa_out;
  2813.      USER_DEFINES *OPTIONS;
  2814. #endif /* HAVE_ANSI */
  2815. {
  2816.   ALLOC_INT index_v;
  2817.   ALLOC_INT index_vv, index_v_vv;
  2818.  
  2819.   fprintf (ptr_asa_out, "\n");
  2820. #if TIME_CALC
  2821.   print_time ("", ptr_asa_out);
  2822. #endif
  2823.  
  2824.   if (OPTIONS->CURVATURE_0 == TRUE)
  2825.     *curvature_flag = FALSE;
  2826.   if (OPTIONS->CURVATURE_0 == -1)
  2827.     *curvature_flag = TRUE;
  2828.  
  2829. #if INT_LONG
  2830.   fprintf (ptr_asa_out,
  2831.       "*index_cost_acceptances = %ld, *current_cost_temperature = %12.7g\n",
  2832.        *index_cost_acceptances, *current_cost_temperature);
  2833.   fprintf (ptr_asa_out,
  2834.        "*accepted_to_generated_ratio = %12.7g,\
  2835.  *number_invalid... = %ld\n",
  2836.        *accepted_to_generated_ratio,
  2837.        (*number_invalid_generated_states));
  2838.   fprintf (ptr_asa_out,
  2839.        "*number_generated = %ld, *number_accepted = %ld\n",
  2840.        *number_generated, *number_accepted);
  2841. #else
  2842.   fprintf (ptr_asa_out,
  2843.        "*index_cost_acceptances = %d, *current_cost_temperature = %12.7g\n",
  2844.        *index_cost_acceptances, *current_cost_temperature);
  2845.   fprintf (ptr_asa_out,
  2846.        "*accepted_to_generated_ratio = %12.7g,\
  2847.  *number_invalid... = %d\n",
  2848.        *accepted_to_generated_ratio,
  2849.        *number_invalid_generated_states);
  2850.   fprintf (ptr_asa_out,
  2851.        "*number_generated = %d, *number_accepted = %d\n",
  2852.        *number_generated, *number_accepted);
  2853. #endif
  2854.  
  2855.   fprintf (ptr_asa_out,
  2856.        "best...->cost = %12.7g,\
  2857.  last...->cost = %12.7g\n",
  2858.        best_generated_state->cost, last_saved_state->cost);
  2859.  
  2860.   /* Note that tangents will not be calculated until reanneal
  2861.      is called, and therefore their listing in the printout only
  2862.      is relevant then */
  2863.  
  2864.   VFOR (index_v)
  2865.   {
  2866.     /* ignore too small ranges */
  2867.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2868.       continue;
  2869.     fprintf (ptr_asa_out,
  2870. #if INT_ALLOC
  2871.          "best_generated_state->parameter[%d] = %12.7g\n\
  2872.  current_user_parameter_temp[%d] = %12.7g\n\
  2873.  tangents[%d]: %12.7g\n",
  2874. #else
  2875. #if INT_LONG
  2876.          "best_generated_state->parameter[%ld] = %12.7g\n\
  2877.  current_user_parameter_temp[%ld] = %12.7g\n\
  2878.  tangents[%ld]: %12.7g\n",
  2879. #else
  2880.          "best_generated_state->parameter[%d] = %12.7g\n\
  2881.  current_user_parameter_temp[%d] = %12.7g\n\
  2882.  tangents[%d]: %12.7g\n",
  2883. #endif
  2884. #endif
  2885.          index_v, best_generated_state->parameter[index_v],
  2886.          index_v, current_user_parameter_temp[index_v],
  2887.          index_v, tangents[index_v]);
  2888.   }
  2889.  
  2890.   if (*curvature_flag == TRUE)
  2891.     {
  2892.       /* print curvatures */
  2893.       VFOR (index_v)
  2894.       {
  2895.     /* ignore too small ranges */
  2896.     if (PARAMETER_RANGE_TOO_SMALL (index_v))
  2897.       continue;
  2898.     VFOR (index_vv)
  2899.     {
  2900.       /* only print upper diagonal of matrix */
  2901.       if (index_v < index_vv)
  2902.         continue;
  2903.       /* ignore too small ranges (index_vv) */
  2904.       if (PARAMETER_RANGE_TOO_SMALL (index_vv))
  2905.         continue;
  2906.  
  2907.       /* index_v_vv: row index_v, column index_vv */
  2908.       index_v_vv = ROW_COL_INDEX (index_v, index_vv);
  2909.  
  2910. #if INT_ALLOC
  2911.       fprintf (ptr_asa_out, "curvature[%d][%d] = %12.7g\n",
  2912. #else
  2913. #if INT_LONG
  2914.       fprintf (ptr_asa_out, "curvature[%ld][%ld] = %12.7g\n",
  2915. #else
  2916.       fprintf (ptr_asa_out, "curvature[%d][%d] = %12.7g\n",
  2917. #endif
  2918. #endif
  2919.            index_v, index_vv, curvature[index_v_vv]);
  2920.     }
  2921.       }
  2922.     }
  2923.   fprintf (ptr_asa_out, "\n");
  2924.   fflush (ptr_asa_out);
  2925.  
  2926. }
  2927.  
  2928. /***********************************************************************
  2929. * print_asa-options
  2930. *    Prints user's selected options
  2931. ***********************************************************************/
  2932. #if HAVE_ANSI
  2933. void
  2934. print_asa_options (FILE * ptr_asa_out, USER_DEFINES * OPTIONS)
  2935. #else
  2936. void
  2937. print_asa_options (ptr_asa_out, OPTIONS)
  2938.      FILE *ptr_asa_out;
  2939.      USER_DEFINES *OPTIONS;
  2940. #endif /* HAVE_ANSI */
  2941. {
  2942.   fprintf (ptr_asa_out,
  2943.        "\t\tADAPTIVE SIMULATED ANNEALING\n\n");
  2944.  
  2945.   fprintf (ptr_asa_out, "%s\n\n", ASA_ID);
  2946.  
  2947.   fprintf (ptr_asa_out, "OPTIONS_FILE = %d\n",
  2948.        (int) OPTIONS_FILE);
  2949.   fprintf (ptr_asa_out, "ASA_LIB = %d\n",
  2950.        (int) ASA_LIB);
  2951.   fprintf (ptr_asa_out, "HAVE_ANSI = %d\n",
  2952.        (int) HAVE_ANSI);
  2953.   fprintf (ptr_asa_out, "IO_PROTOTYPES = %d\n",
  2954.        (int) IO_PROTOTYPES);
  2955.   fprintf (ptr_asa_out, "TIME_CALC = %d\n",
  2956.        (int) TIME_CALC);
  2957.   fprintf (ptr_asa_out, "TIME_STD = %d\n",
  2958.        (int) TIME_STD);
  2959.   fprintf (ptr_asa_out, "INT_LONG = %d\n",
  2960.        (int) INT_LONG);
  2961.   fprintf (ptr_asa_out, "INT_ALLOC = %d\n",
  2962.        (int) INT_ALLOC);
  2963.   fprintf (ptr_asa_out, "SMALL_FLOAT = %g\n",
  2964.        (double) SMALL_FLOAT);
  2965.   fprintf (ptr_asa_out, "MIN_DOUBLE = %g\n",
  2966.        (double) MIN_DOUBLE);
  2967.   fprintf (ptr_asa_out, "MAX_DOUBLE = %g\n",
  2968.        (double) MAX_DOUBLE);
  2969.   fprintf (ptr_asa_out, "EPS_DOUBLE = %g\n",
  2970.        (double) EPS_DOUBLE);
  2971.   fprintf (ptr_asa_out, "NO_PARAM_TEMP_TEST = %d\n",
  2972.        (int) NO_PARAM_TEMP_TEST);
  2973.   fprintf (ptr_asa_out, "NO_COST_TEMP_TEST = %d\n",
  2974.        (int) NO_COST_TEMP_TEST);
  2975.   fprintf (ptr_asa_out, "SELF_OPTIMIZE = %d\n",
  2976.        (int) SELF_OPTIMIZE);
  2977.   fprintf (ptr_asa_out, "ASA_TEST = %d\n",
  2978.        (int) ASA_TEST);
  2979.   fprintf (ptr_asa_out, "ASA_TEMPLATE = %d\n",
  2980.        (int) ASA_TEMPLATE);
  2981.   fprintf (ptr_asa_out, "OPTIONAL_DATA = %d\n",
  2982.        (int) OPTIONAL_DATA);
  2983.   fprintf (ptr_asa_out, "USER_COST_SCHEDULE = %d\n",
  2984.        (int) USER_COST_SCHEDULE);
  2985.   fprintf (ptr_asa_out, "USER_REANNEAL_FUNCTION = %d\n",
  2986.        (int) USER_REANNEAL_FUNCTION);
  2987.   fprintf (ptr_asa_out, "ASA_SAMPLE = %d\n",
  2988.        (int) ASA_SAMPLE);
  2989.   fprintf (ptr_asa_out, "ASA_PARALLEL = %d\n\n",
  2990.        (int) ASA_PARALLEL);
  2991.  
  2992.   fprintf (ptr_asa_out, "ASA_PRINT = %d\n",
  2993.        (int) ASA_PRINT);
  2994.   fprintf (ptr_asa_out, "ASA_OUT = %s\n",
  2995. #if USER_ASA_OUT
  2996.        OPTIONS->asa_out_file);
  2997. #else
  2998.        ASA_OUT);
  2999. #endif
  3000.   fprintf (ptr_asa_out, "USER_ASA_OUT = %d\n",
  3001.        (int) USER_ASA_OUT);
  3002.   fprintf (ptr_asa_out, "ASA_PRINT_INTERMED = %d\n",
  3003.        (int) ASA_PRINT_INTERMED);
  3004.   fprintf (ptr_asa_out, "ASA_PRINT_MORE = %d\n\n",
  3005.        (int) ASA_PRINT_MORE);
  3006.  
  3007. #if INT_LONG
  3008.   fprintf (ptr_asa_out, "OPTIONS->LIMIT_ACCEPTANCES = %ld\n",
  3009.        (LONG_INT) OPTIONS->LIMIT_ACCEPTANCES);
  3010.   fprintf (ptr_asa_out, "OPTIONS->LIMIT_GENERATED = %ld\n",
  3011.        (LONG_INT) OPTIONS->LIMIT_GENERATED);
  3012. #else
  3013.   fprintf (ptr_asa_out, "OPTIONS->LIMIT_ACCEPTANCES = %d\n",
  3014.        (LONG_INT) OPTIONS->LIMIT_ACCEPTANCES);
  3015.   fprintf (ptr_asa_out, "OPTIONS->LIMIT_GENERATED = %d\n",
  3016.        (LONG_INT) OPTIONS->LIMIT_GENERATED);
  3017. #endif
  3018.   fprintf (ptr_asa_out, "OPTIONS->LIMIT_INVALID_GENERATED_STATES = %d\n",
  3019.        OPTIONS->LIMIT_INVALID_GENERATED_STATES);
  3020.   fprintf (ptr_asa_out, "OPTIONS->ACCEPTED_TO_GENERATED_RATIO = %g\n\n",
  3021.        OPTIONS->ACCEPTED_TO_GENERATED_RATIO);
  3022.  
  3023.   fprintf (ptr_asa_out, "OPTIONS->COST_PRECISION = %g\n",
  3024.        OPTIONS->COST_PRECISION);
  3025.   fprintf (ptr_asa_out, "OPTIONS->MAXIMUM_COST_REPEAT = %d\n",
  3026.        OPTIONS->MAXIMUM_COST_REPEAT);
  3027.   fprintf (ptr_asa_out, "OPTIONS->NUMBER_COST_SAMPLES = %d\n",
  3028.        OPTIONS->NUMBER_COST_SAMPLES);
  3029.   fprintf (ptr_asa_out, "OPTIONS->TEMPERATURE_RATIO_SCALE = %g\n",
  3030.        OPTIONS->TEMPERATURE_RATIO_SCALE);
  3031.   fprintf (ptr_asa_out, "OPTIONS->COST_PARAMETER_SCALE = %g\n",
  3032.        OPTIONS->COST_PARAMETER_SCALE);
  3033.   fprintf (ptr_asa_out, "OPTIONS->TEMPERATURE_ANNEAL_SCALE = %g\n",
  3034.        OPTIONS->TEMPERATURE_ANNEAL_SCALE);
  3035.   fprintf (ptr_asa_out, "OPTIONS->USER_INITIAL_COST_TEMP = %d\n\n",
  3036.        OPTIONS->USER_INITIAL_COST_TEMP);
  3037.  
  3038.   fprintf (ptr_asa_out, "OPTIONS->INCLUDE_INTEGER_PARAMETERS = %d\n",
  3039.        OPTIONS->INCLUDE_INTEGER_PARAMETERS);
  3040.   fprintf (ptr_asa_out, "OPTIONS->USER_INITIAL_PARAMETERS = %d\n",
  3041.        OPTIONS->USER_INITIAL_PARAMETERS);
  3042. #if INT_ALLOC
  3043.   fprintf (ptr_asa_out, "OPTIONS->SEQUENTIAL_PARAMETERS = %d\n",
  3044.        (int) OPTIONS->SEQUENTIAL_PARAMETERS);
  3045. #else
  3046. #if INT_LONG
  3047.   fprintf (ptr_asa_out, "OPTIONS->SEQUENTIAL_PARAMETERS = %ld\n",
  3048.        (LONG_INT) OPTIONS->SEQUENTIAL_PARAMETERS);
  3049. #else
  3050.   fprintf (ptr_asa_out, "OPTIONS->SEQUENTIAL_PARAMETERS = %d\n",
  3051.        (LONG_INT) OPTIONS->SEQUENTIAL_PARAMETERS);
  3052. #endif
  3053. #endif
  3054.   fprintf (ptr_asa_out, "OPTIONS->INITIAL_PARAMETER_TEMPERATURE = %g\n",
  3055.        OPTIONS->INITIAL_PARAMETER_TEMPERATURE);
  3056.   fprintf (ptr_asa_out, "OPTIONS->RATIO_TEMPERATURE_SCALES = %d\n",
  3057.        OPTIONS->RATIO_TEMPERATURE_SCALES);
  3058.   fprintf (ptr_asa_out, "OPTIONS->USER_INITIAL_PARAMETERS_TEMPS = %d\n\n",
  3059.        OPTIONS->USER_INITIAL_PARAMETERS_TEMPS);
  3060.  
  3061.   fprintf (ptr_asa_out, "OPTIONS->TESTING_FREQUENCY_MODULUS = %d\n",
  3062.        OPTIONS->TESTING_FREQUENCY_MODULUS);
  3063.   fprintf (ptr_asa_out, "OPTIONS->ACTIVATE_REANNEAL = %d\n",
  3064.        OPTIONS->ACTIVATE_REANNEAL);
  3065.   fprintf (ptr_asa_out, "OPTIONS->REANNEAL_RESCALE = %g\n",
  3066.        OPTIONS->REANNEAL_RESCALE);
  3067. #if INT_LONG
  3068.   fprintf (ptr_asa_out, "OPTIONS->MAXIMUM_REANNEAL_INDEX = %ld\n",
  3069.        (LONG_INT) OPTIONS->MAXIMUM_REANNEAL_INDEX);
  3070. #else
  3071.   fprintf (ptr_asa_out, "OPTIONS->MAXIMUM_REANNEAL_INDEX = %d\n\n",
  3072.        OPTIONS->MAXIMUM_REANNEAL_INDEX);
  3073. #endif
  3074.  
  3075.   fprintf (ptr_asa_out, "OPTIONS->DELTA_X = %g\n",
  3076.        OPTIONS->DELTA_X);
  3077.   fprintf (ptr_asa_out, "OPTIONS->DELTA_PARAMETERS = %d\n",
  3078.        OPTIONS->DELTA_PARAMETERS);
  3079.   fprintf (ptr_asa_out, "OPTIONS->USER_TANGENTS = %d\n",
  3080.        OPTIONS->USER_TANGENTS);
  3081.   fprintf (ptr_asa_out, "OPTIONS->CURVATURE_0 = %d\n\n",
  3082.        OPTIONS->CURVATURE_0);
  3083.  
  3084.   fprintf (ptr_asa_out, "OPTIONS->QUENCH_PARAMETERS = %d\n",
  3085.        OPTIONS->QUENCH_PARAMETERS);
  3086.   fprintf (ptr_asa_out, "OPTIONS->QUENCH_COST = %d\n\n",
  3087.        OPTIONS->QUENCH_COST);
  3088.  
  3089.   fprintf (ptr_asa_out, "\n");
  3090. }
  3091.  
  3092. #if TIME_CALC
  3093. /***********************************************************************
  3094. * print_time
  3095. *    This calculates the time and runtime and prints it.
  3096. ***********************************************************************/
  3097. #if HAVE_ANSI
  3098. void
  3099. print_time (char *message, FILE * ptr_asa_out)
  3100. #else
  3101. void
  3102. print_time (message, ptr_asa_out)
  3103.      char *message;
  3104.      FILE *ptr_asa_out;
  3105. #endif /* HAVE_ANSI */
  3106. {
  3107.   int who = RUSAGE_SELF;    /* Check our own time */
  3108.   struct rusage usage;
  3109.  
  3110.   /* get the resource usage information */
  3111. #if TIME_STD
  3112.   syscall (SYS_GETRUSAGE, who, &usage);
  3113. #else
  3114.   getrusage (who, &usage);
  3115. #endif
  3116.  
  3117.   /* print the usage time in reasonable form */
  3118.   aux_print_time (&usage.ru_utime, message, ptr_asa_out);
  3119. }
  3120.  
  3121. /***********************************************************************
  3122. * aux_print_time
  3123. *      auxiliary print the time routine
  3124. ***********************************************************************/
  3125. #if HAVE_ANSI
  3126. void
  3127. aux_print_time (struct timeval *time, char *message,
  3128.         FILE * ptr_asa_out)
  3129. #else
  3130. void
  3131. aux_print_time (time, message, ptr_asa_out)
  3132.      struct timeval *time;
  3133.      char *message;
  3134.      FILE *ptr_asa_out;
  3135. #endif /* HAVE_ANSI */
  3136. {
  3137.   static double sx;
  3138.   double us, s, m, h;
  3139.   double ds, dm, dh;
  3140.  
  3141.   /* calculate the new microseconds, seconds, minutes, hours
  3142.      and the differences since the last call */
  3143.   us = (double) ((int) ((double) EPS_DOUBLE + time->tv_usec)) / 1.E6;
  3144.   s = (double) ((int) ((double) EPS_DOUBLE + time->tv_sec)) + us;
  3145.   ds = s - sx;
  3146.   sx = s;
  3147.  
  3148.   h = (int) ((double) EPS_DOUBLE + s / 3600.);
  3149.   m = (int) ((double) EPS_DOUBLE + s / 60.) - 60. * h;
  3150.   s -= (3600. * h + 60. * m);
  3151.   dh = (int) ((double) EPS_DOUBLE + ds / 3600.);
  3152.   dm = (int) ((double) EPS_DOUBLE + ds / 60.) - 60. * dh;
  3153.   ds -= (3600. * dh + 60. * dm);
  3154.  
  3155.   /* print the statistics */
  3156.   fprintf (ptr_asa_out,
  3157.        "%s:time: %gh %gm %gs; incr: %gh %gm %gs\n",
  3158.        message, h, m, s, dh, dm, ds);
  3159. }
  3160.  
  3161. #endif /* TIME_CALC */
  3162. #endif /* ASA_PRINT */
  3163.